Matrix multiplication with cblas, block works now, reformat lapack cheat a bit

Author: cirosantilli

lapack/README.md

@@ -8,7 +8,7 @@ BLAS and LAPACK are:
 
 -   de-facto standards
 
--   non-parallel 
+-   non-parallel
 
 -   originally written in Fortran
 
@@ -30,8 +30,6 @@ BLAS contains low level functions such as:
 - vector matrix multiplication
 - matrix matrix multiplication
 
-The BLAS project provides `cblas.h`, which contains a C interface for BLAS.
-
 ### LAPACK
 
 LAPACK contains higher level functions such as:
@@ -61,21 +59,19 @@ Implements full BLAS, but only part of LAPACK.
 
 Has C interface.
 
-## Installation on Ubuntu
-
-### Fortran
-
-    sudo aptitude install liblapack-dev liblapack-doc libblas-doc
+## C interface
 
-### C interface
+The BLAS project provides `cblas.h`, which contains a C interface for BLAS (TODO but also an implementation?)
 
-via atlas:
+Via atlas:
 
-    sudo aptitude install 
+    sudo aptitude install libatlas-dev
+    gcc -lcblas
 
-via LAPACKE (`libblas-dev` already contains `cblas.h`):
+Via LAPACKE (`libblas-dev` already contains `cblas.h`):
 
     sudo aptitude install liblapacke-dev
+    gcc -lblas
 
 ## Levels
 
@@ -85,6 +81,8 @@ via LAPACKE (`libblas-dev` already contains `cblas.h`):
 
 ## Function naming conventions
 
+<http://www.netlib.org/lapack/lug/node24.html>
+
 The functions are named according to the pattern:
 
     XYYZZZ
@@ -108,6 +106,8 @@ Where:
 -   `ZZ`: computation to be done:
 
     - `SV`: SolVe linear system
+    - `MM`: Matrix Multiply
+    - `LS`: Least Squares (overdetermined system)
 
 ## Sources
 
@@ -118,3 +118,7 @@ Where:
 - user's guide. algorithm info <http://www.netlib.org/lapack/lug/>
 
 - <http://www.tat.physik.uni-tuebingen.de/~kley/lehre/ftn77/tutorial/blas.html>
+
+## LAPACKE
+
+<http://stackoverflow.com/questions/26875415/difference-between-lapacke-and-lapack>

lapack/c.c

@@ -14,103 +14,97 @@ interfaces provided by their respective projects, repectively through
 #include <lapacke.h>
 
 /**
-* assert two integers are equal
-* if not, print them to stderr and assert false
-*/
-void assert_eqi( int i1, int i2 )
-{
-    if( i1 != i2 )
-    {
-        fprintf( stderr, "%d\n%d\n\n\n", i1, i2 );
-        assert( false );
+ * assert two integers are equal
+ * if not, print them to stderr and assert false
+ */
+void assert_eqi(int i1, int i2) {
+    if (i1 != i2) {
+        fprintf(stderr, "%d\n%d\n\n\n", i1, i2);
+        assert(false);
     }
 }
 
-/*
-* assert two doubles are equal within err precision
-* if not, print them to stderr
-*/
-void assert_eqd( double d1, double d2, double err )
-{
-    if( fabs( d1 - d2 ) > err )
-    {
-        fprintf( stderr, "%f\n%f\n\n\n", d1, d2 );
-        assert( false );
+/**
+ * assert two doubles are equal within err precision
+ * if not, print them to stderr
+ */
+void assert_eqd(double d1, double d2, double err) {
+    if (fabs(d1 - d2) > err) {
+        fprintf(stderr, "%f\n%f\n\n\n", d1, d2);
+        assert(false);
     }
 }
 
-/*
-* print an array of doubles to stderr
-*/
-void print_vecd( int n, double * v )
-{
+/** print an array of doubles to stderr */
+void print_vecd(int n, double * v) {
     int i;
-    for(i=0; i<n; i++)
-    {
-        fprintf( stderr, "%.2f ", v[i] );
+    for (i=0; i<n; i++) {
+        fprintf(stderr, "%.2f ", v[i]);
     }
     printf("\n");
 }
 
 /*
-bool eq_vecd2( int n, double * v1, double * v2, double err)
+bool eq_vecd2(int n, double * v1, double * v2, double err)
 {
-    if( fabs( d1 - d2 ) > err )
+    if(fabs(d1 - d2) > err)
         return false
     return true;
 }
 
-void assert_eqvd( int n, double * v1, double * v2, double err){
-    if( eq_vecd( n, v1, v2, err ) != true ){
-        print_vecd( n, v1 );
-        print_vecd( n, v2 );
+void assert_eqvd(int n, double * v1, double * v2, double err){
+    if(eq_vecd(n, v1, v2, err) != true){
+        print_vecd(n, v1);
+        print_vecd(n, v2);
     }
 }
 */
 
-int main(void)
-{
+int main(void) {
     int info, ipiv2[2];
     float err = 1e-6;
     float x2[2], b2[2], c2[2];
     float a2x2[2][2];
 
-    //#cblas
-
+	/* cblas */
+	{
         x2[0] =  1.0;
         x2[1] = -2.0;
-        assert_eqd( cblas_snrm2(2, x2, 1), sqrt(5.0), err );
-
-    //#lapacke
-
-        //1 2 x = 5
-        //3 4 y   11
-
-        //x = 1
-        //y = 1
-
-        a2x2[0][0] = 1.0;
-        a2x2[1][0] = 2.0;
-        a2x2[0][1] = 3.0;
-        a2x2[1][1] = 4.0;
-        b2[0] = 5.;
-        b2[1] = 11.;
+        assert_eqd(cblas_snrm2(2, x2, 1), sqrt(5.0), err);
+    }
 
-        info = LAPACKE_sgesv(
-            LAPACK_COL_MAJOR, //COL or ROW
-            2,                //n
-            1,                //nrhs
-            &a2x2[0][0],
-            2,                //lda
-            &ipiv2[0],
-            &b2[0],
-            2                 //ldb
-        );
-        c2[0] = 1.0;
-        c2[1] = 2.0;
-        assert_eqi( info, 0 );
-        assert_eqd( b2[0], c2[0], err );
-        assert_eqd( b2[1], c2[1], err );
+	/* lapacke */
+	{
+		/*
+		sgesv
+
+		Matrix vector multiply.
+		*/
+		{
+			a2x2[0][0] = 1.0;
+			a2x2[1][0] = 2.0;
+			a2x2[0][1] = 3.0;
+			a2x2[1][1] = 4.0;
+			b2[0] = 5.;
+			b2[1] = 11.;
+
+			info = LAPACKE_sgesv(
+				LAPACK_COL_MAJOR,
+				2,
+				1,
+				&a2x2[0][0],
+				2,
+				&ipiv2[0],
+				&b2[0],
+				2
+			);
+			c2[0] = 1.0;
+			c2[1] = 2.0;
+			assert_eqi(info, 0);
+			assert_eqd(b2[0], c2[0], err);
+			assert_eqd(b2[1], c2[1], err);
+		}
+    }
 
     return EXIT_SUCCESS;
 }

lapack/configure

@@ -1,7 +1,7 @@
 #!/usr/bin/env bash
 sudo aptitude update
 # BLAS C/Fotran and LAPACK Fortran:
-sudo aptitude install -y liblapack-dev
+sudo aptitude install -y liblapack-dev liblapacke-dev
 # Lapack C via LAPACKE:
 # TODO: how to install on Ubuntu 12.04?
 # The following works only on later Ubuntu

lapack/main.f

@@ -171,7 +171,7 @@ program main
             !return status returned on `info`:
             !pivots returned on `pivots`:
             call sgesv( n, nrhs, a2x2, lda, pivots, b2, ldb, info )
-            
+
             c2(1) = 1.0
             c2(2) = 2.0
             call assert_eqi( info, 0 )

opencl/Makefile_params

@@ -1 +1,3 @@
-LIBS := -lm -pthread -lOpenCL
+LIBS := -lm -pthread -lOpenCL -lblas
+# blas from libatlas-base-dev
+# cblas from libblas-dev