Mercurial > octave-nkf

changeset 5983:ae09df27153f

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2006-09-12 02:15:47 by jwe]
author jwe
date Tue, 12 Sep 2006 02:15:47 +0000
parents 1ed991f0ed61
children 82a73f5dadd9
files libcruft/ChangeLog libcruft/blas-xtra/xddot.f libcruft/blas-xtra/xzdotu.f liboctave/CMatrix.cc liboctave/CRowVector.cc liboctave/ChangeLog liboctave/chNDArray.cc liboctave/dMatrix.cc liboctave/dRowVector.cc
diffstat 9 files changed, 107 insertions(+), 25 deletions(-) [+]
line wrap: on
line diff
--- a/libcruft/ChangeLog	Tue Sep 05 20:13:50 2006 +0000
+++ b/libcruft/ChangeLog	Tue Sep 12 02:15:47 2006 +0000
@@ -1,3 +1,7 @@
+2006-09-11  John W. Eaton  <jwe@octave.org>
+
+	* blas-xtra/xddot.f, blas-xtra/xzdotu.f: New files.
+
 2006-06-01  David Bateman  <dbateman@free.fr>
 
 	* slatec-fn/dpchim.f, slatec-fn/dpchst.f: New files.
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/libcruft/blas-xtra/xddot.f	Tue Sep 12 02:15:47 2006 +0000
@@ -0,0 +1,6 @@
+      subroutine xddot (n, dx, incx, dy, incy, retval)
+      double precision ddot, dx(*), dy(*), retval
+      integer incx, incy
+      retval = ddot (n, dx, incx, dy, incy)
+      return
+      end
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/libcruft/blas-xtra/xzdotu.f	Tue Sep 12 02:15:47 2006 +0000
@@ -0,0 +1,6 @@
+      subroutine xzdotu (n, zx, incx, zy, incy, retval)
+      double complex zdotu, zx(*), zy(*), retval
+      integer incx, incy
+      retval = zdotu (n, zx, incx, zy, incy)
+      return
+      end
--- a/liboctave/CMatrix.cc	Tue Sep 05 20:13:50 2006 +0000
+++ b/liboctave/CMatrix.cc	Tue Sep 12 02:15:47 2006 +0000
@@ -86,6 +86,17 @@
 			   F77_CHAR_ARG_LEN_DECL);
 
   F77_RET_T
+  F77_FUNC (zgemv, ZGEMV) (F77_CONST_CHAR_ARG_DECL,
+                           const octave_idx_type&, const octave_idx_type&, const Complex&,
+                           const Complex*, const octave_idx_type&, const Complex*,
+                           const octave_idx_type&, const Complex&, Complex*, const octave_idx_type&
+                           F77_CHAR_ARG_LEN_DECL);
+
+  F77_RET_T
+  F77_FUNC (xzdotu, XZDOTU) (const octave_idx_type&, const Complex*, const octave_idx_type&,
+			     const Complex*, const octave_idx_type&, Complex&);
+
+  F77_RET_T
   F77_FUNC (zgetrf, ZGETRF) (const octave_idx_type&, const octave_idx_type&, Complex*, const octave_idx_type&,
 			     octave_idx_type*, octave_idx_type&);
 
@@ -3619,12 +3630,23 @@
 	  retval.resize (nr, a_nc);
 	  Complex *c = retval.fortran_vec ();
 
-	  F77_XFCN (zgemm, ZGEMM, (F77_CONST_CHAR_ARG2 ("N", 1),
-				   F77_CONST_CHAR_ARG2 ("N", 1),
-				   nr, a_nc, nc, 1.0, m.data (),
-				   ld, a.data (), lda, 0.0, c, nr
-				   F77_CHAR_ARG_LEN (1)
-				   F77_CHAR_ARG_LEN (1)));
+	  if (a_nc == 1)
+	    {
+	      if (nr == 1)
+		F77_FUNC (xzdotu, XZDOTU) (nc, m.data (), 1, a.data (), 1, *c);
+	      else
+		F77_XFCN (zgemv, ZGEMV, (F77_CONST_CHAR_ARG2 ("N", 1),
+					 nr, nc, 1.0,  m.data (), ld,
+					 a.data (), 1, 0.0, c, 1
+					 F77_CHAR_ARG_LEN (1)));
+	    }
+	  else
+	    F77_XFCN (zgemm, ZGEMM, (F77_CONST_CHAR_ARG2 ("N", 1),
+				     F77_CONST_CHAR_ARG2 ("N", 1),
+				     nr, a_nc, nc, 1.0, m.data (),
+				     ld, a.data (), lda, 0.0, c, nr
+				     F77_CHAR_ARG_LEN (1)
+				     F77_CHAR_ARG_LEN (1)));
 
 	  if (f77_exception_encountered)
 	    (*current_liboctave_error_handler)
--- a/liboctave/CRowVector.cc	Tue Sep 05 20:13:50 2006 +0000
+++ b/liboctave/CRowVector.cc	Tue Sep 12 02:15:47 2006 +0000
@@ -45,6 +45,10 @@
 			   const Complex*, const octave_idx_type&, const Complex*,
 			   const octave_idx_type&, const Complex&, Complex*, const octave_idx_type&
 			   F77_CHAR_ARG_LEN_DECL);
+
+  F77_RET_T
+  F77_FUNC (xzdotu, XZDOTU) (const octave_idx_type&, const Complex*, const octave_idx_type&,
+			     const Complex*, const octave_idx_type&, Complex&);
 }
 
 // Complex Row Vector class
@@ -488,18 +492,16 @@
 Complex
 operator * (const ComplexRowVector& v, const ComplexColumnVector& a)
 {
+  Complex retval (0.0, 0.0);
+
   octave_idx_type len = v.length ();
 
   octave_idx_type a_len = a.length ();
 
   if (len != a_len)
-    {
-      gripe_nonconformant ("operator *", len, a_len);
-      return 0.0;
-    }
-
-  Complex retval (0.0, 0.0);
-
+    gripe_nonconformant ("operator *", len, a_len);
+  else if (len != 0)
+    F77_FUNC (xzdotu, XZDOTU) (len, v.data (), 1, a.data (), 1, retval);
   for (octave_idx_type i = 0; i < len; i++)
     retval += v.elem (i) * a.elem (i);
 
--- a/liboctave/ChangeLog	Tue Sep 05 20:13:50 2006 +0000
+++ b/liboctave/ChangeLog	Tue Sep 12 02:15:47 2006 +0000
@@ -1,3 +1,22 @@
+2006-09-11  John W. Eaton  <jwe@octave.org>
+
+	* dMatrix.cc (operator * (const Matrix&, const Matrix&))):
+	Handle M*v and rv*cv special cases. 
+	* CMatrix.cc (operator * (const ComplexMatrix&, const
+	ComplexMatrix&))): Likewise.
+	From Luis F. Ortiz <lortiz@interactivesupercomputing.com>.
+
+	* dRowVector.cc (operator * (const RowVector&, const
+	ColumnVector&)): Call xddot here instead of using a Fortran
+	function directly.
+	* CRowVector.cc (operator * (const ComplexRowVector&, const
+	ComplexColumnVector&)): Call xzdotu here.
+
+2006-09-05  John W. Eaton  <jwe@octave.org>
+
+	* chNDArray.cc (charNDArray::any, charNDArray::all): Compare
+	elements to '\0', not ' '.
+
 2006-08-25  John W. Eaton  <jwe@octave.org>
 
 	* mx-inlines.cc (MX_ND_REDUCTION): Special case for 0x0 arrays.
--- a/liboctave/chNDArray.cc	Tue Sep 05 20:13:50 2006 +0000
+++ b/liboctave/chNDArray.cc	Tue Sep 12 02:15:47 2006 +0000
@@ -37,13 +37,13 @@
 boolNDArray
 charNDArray::all (int dim) const
 {
-  MX_ND_ANY_ALL_REDUCTION (MX_ND_ALL_EVAL (elem (iter_idx) == ' '), true);
+  MX_ND_ANY_ALL_REDUCTION (MX_ND_ALL_EVAL (elem (iter_idx) == '\0'), true);
 }
 
 boolNDArray
 charNDArray::any (int dim) const
 {
-  MX_ND_ANY_ALL_REDUCTION (MX_ND_ANY_EVAL (elem (iter_idx) != ' '), false);
+  MX_ND_ANY_ALL_REDUCTION (MX_ND_ANY_EVAL (elem (iter_idx) != '\0'), false);
 }
 
 charNDArray
--- a/liboctave/dMatrix.cc	Tue Sep 05 20:13:50 2006 +0000
+++ b/liboctave/dMatrix.cc	Tue Sep 12 02:15:47 2006 +0000
@@ -83,6 +83,18 @@
 			   F77_CHAR_ARG_LEN_DECL);
 
   F77_RET_T
+  F77_FUNC (dgemv, DGEMV) (F77_CONST_CHAR_ARG_DECL,
+			   const octave_idx_type&, const octave_idx_type&, const double&,
+			   const double*, const octave_idx_type&, const double*,
+			   const octave_idx_type&, const double&, double*,
+			   const octave_idx_type&
+			   F77_CHAR_ARG_LEN_DECL);
+
+  F77_RET_T
+  F77_FUNC (xddot, XDDOT) (const octave_idx_type&, const double*, const octave_idx_type&,
+			   const double*, const octave_idx_type&, double&);
+
+  F77_RET_T
   F77_FUNC (dgetrf, DGETRF) (const octave_idx_type&, const octave_idx_type&, double*, const octave_idx_type&,
 		      octave_idx_type*, octave_idx_type&);
 
@@ -2993,12 +3005,23 @@
 	  retval.resize (nr, a_nc);
 	  double *c = retval.fortran_vec ();
 
-	  F77_XFCN (dgemm, DGEMM, (F77_CONST_CHAR_ARG2 ("N", 1),
-				   F77_CONST_CHAR_ARG2 ("N", 1),
-				   nr, a_nc, nc, 1.0, m.data (),
-				   ld, a.data (), lda, 0.0, c, nr
-				   F77_CHAR_ARG_LEN (1)
-				   F77_CHAR_ARG_LEN (1)));
+	  if (a_nc == 1)
+	    {
+	      if (nr == 1)
+		F77_FUNC (xddot, XDDOT) (nc, m.data (), 1, a.data (), 1, *c);
+	      else
+		F77_XFCN (dgemv, DGEMV, (F77_CONST_CHAR_ARG2 ("N", 1),
+					 nr, nc, 1.0,  m.data (), ld,
+					 a.data (), 1, 0.0, c, 1
+					 F77_CHAR_ARG_LEN (1)));
+            }
+	  else
+	    F77_XFCN (dgemm, DGEMM, (F77_CONST_CHAR_ARG2 ("N", 1),
+				     F77_CONST_CHAR_ARG2 ("N", 1),
+				     nr, a_nc, nc, 1.0, m.data (),
+				     ld, a.data (), lda, 0.0, c, nr
+				     F77_CHAR_ARG_LEN (1)
+				     F77_CHAR_ARG_LEN (1)));
 
 	  if (f77_exception_encountered)
 	    (*current_liboctave_error_handler)
--- a/liboctave/dRowVector.cc	Tue Sep 05 20:13:50 2006 +0000
+++ b/liboctave/dRowVector.cc	Tue Sep 12 02:15:47 2006 +0000
@@ -45,9 +45,9 @@
 			   const double*, const octave_idx_type&, const double*,
 			   const octave_idx_type&, const double&, double*, const octave_idx_type&
 			   F77_CHAR_ARG_LEN_DECL);
-
-  double F77_FUNC (ddot, DDOT) (const octave_idx_type&, const double*, const octave_idx_type&,
-				const double*, const octave_idx_type&);
+  F77_RET_T
+  F77_FUNC (xddot, XDDOT) (const octave_idx_type&, const double*, const octave_idx_type&,
+			   const double*, const octave_idx_type&, double&);
 }
 
 // Row Vector class.
@@ -365,7 +365,7 @@
   if (len != a_len)
     gripe_nonconformant ("operator *", len, a_len);
   else if (len != 0)
-    retval = F77_FUNC (ddot, DDOT) (len, v.data (), 1, a.data (), 1);
+    F77_FUNC (xddot, XDDOT) (len, v.data (), 1, a.data (), 1, retval);
 
   return retval;
 }