--- Numeric-24.2/Packages/dotblas/dotblas/_dotblas.c.orig	2005-04-07 06:15:47.000000000 +0100
+++ Numeric-24.2/Packages/dotblas/dotblas/_dotblas.c	2007-05-25 17:39:30.382025500 +0100
@@ -52,7 +52,7 @@
 static PyObject *dotblas_matrixproduct(PyObject *dummy, PyObject *args) {
   PyObject *op1, *op2;
   PyArrayObject *ap1, *ap2, *ret;
-  int i, j, l, lda, ldb, matchDim = -1, otherDim = -1;
+  int i, j, l, lda, ldb, ldc, matchDim = -1, otherDim = -1;
   int typenum;
   int dimensions[MAX_DIMS], nd;
   static const float oneF[2] = {1.0, 0.0};
@@ -248,33 +248,34 @@
     /* Matrix matrix multiplication -- Level 3 BLAS */  
     lda = (ap1->dimensions[1] > 1 ? ap1->dimensions[1] : 1);
     ldb = (ap2->dimensions[1] > 1 ? ap2->dimensions[1] : 1);
+    ldc = (ret->dimensions[1] > 1 ? ret->dimensions[1] : 1);
     if (typenum == PyArray_DOUBLE) {
       cblas_dgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans,
 		  ap1->dimensions[0], ap2->dimensions[1], ap2->dimensions[0],
 		  1.0, (double *)ap1->data, lda,
 		  (double *)ap2->data, ldb,
-		  0.0, (double *)ret->data, ldb);
+		  0.0, (double *)ret->data, ldc);
     }
     else if (typenum == PyArray_FLOAT) {
       cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans,
 		  ap1->dimensions[0], ap2->dimensions[1], ap2->dimensions[0],
 		  1.0, (float *)ap1->data, lda,
 		  (float *)ap2->data, ldb,
-		  0.0, (float *)ret->data, ldb);
+		  0.0, (float *)ret->data, ldc);
     }
     else if (typenum == PyArray_CDOUBLE) {
       cblas_zgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans,
 		  ap1->dimensions[0], ap2->dimensions[1], ap2->dimensions[0],
 		  oneD, (double *)ap1->data, lda,
 		  (double *)ap2->data, ldb,
-		  zeroD, (double *)ret->data, ldb);
+		  zeroD, (double *)ret->data, ldc);
     }
     else if (typenum == PyArray_CFLOAT) {
       cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans,
 		  ap1->dimensions[0], ap2->dimensions[1], ap2->dimensions[0],
 		  oneF, (float *)ap1->data, lda,
 		  (float *)ap2->data, ldb,
-		  zeroF, (float *)ret->data, ldb);
+		  zeroF, (float *)ret->data, ldc);
     }
   }
   else {
@@ -334,7 +335,7 @@
 static PyObject *dotblas_innerproduct(PyObject *dummy, PyObject *args) {
   PyObject *op1, *op2;
   PyArrayObject *ap1, *ap2, *ret;
-  int i, j, l, lda, ldb;
+  int i, j, l, lda, ldb, ldc;
   int typenum;
   int dimensions[MAX_DIMS], nd;
   static const float oneF[2] = {1.0, 0.0};
@@ -525,33 +526,34 @@
     /* Matrix matrix multiplication -- Level 3 BLAS */  
     lda = (ap1->dimensions[1] > 1 ? ap1->dimensions[1] : 1);
     ldb = (ap2->dimensions[1] > 1 ? ap2->dimensions[1] : 1);
+    ldc = (ret->dimensions[1] > 1 ? ret->dimensions[1] : 1);
     if (typenum == PyArray_DOUBLE) {
       cblas_dgemm(CblasRowMajor, CblasNoTrans, CblasTrans,
 		  ap1->dimensions[0], ap2->dimensions[0], ap1->dimensions[1],
 		  1.0, (double *)ap1->data, lda,
 		  (double *)ap2->data, ldb,
-		  0.0, (double *)ret->data, ldb);
+		  0.0, (double *)ret->data, ldc);
     }
     else if (typenum == PyArray_FLOAT) {
       cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasTrans,
 		  ap1->dimensions[0], ap2->dimensions[0], ap1->dimensions[1],
 		  1.0, (float *)ap1->data, lda,
 		  (float *)ap2->data, ldb,
-		  0.0, (float *)ret->data, ldb);
+		  0.0, (float *)ret->data, ldc);
     }
     else if (typenum == PyArray_CDOUBLE) {
       cblas_zgemm(CblasRowMajor, CblasNoTrans, CblasTrans,
 		  ap1->dimensions[0], ap2->dimensions[0], ap1->dimensions[1],
 		  oneD, (double *)ap1->data, lda,
 		  (double *)ap2->data, ldb,
-		  zeroD, (double *)ret->data, ldb);
+		  zeroD, (double *)ret->data, ldc);
     }
     else if (typenum == PyArray_CFLOAT) {
       cblas_cgemm(CblasRowMajor, CblasNoTrans, CblasTrans,
 		  ap1->dimensions[0], ap2->dimensions[0], ap1->dimensions[1],
 		  oneF, (float *)ap1->data, lda,
 		  (float *)ap2->data, ldb,
-		  zeroF, (float *)ret->data, ldb);
+		  zeroF, (float *)ret->data, ldc);
     }
   }
   else {