From 324e0a33dd7a3bf5f47c88b37de68dfca917ef2d Mon Sep 17 00:00:00 2001
From: Joseph Redmon <pjreddie@gmail.com>
Date: Fri, 18 Apr 2014 06:14:13 +0000
Subject: [PATCH] Better alternating between video cards

---
 src/mini_blas.c |  201 ++++++++++++++++++++++++++++++++++++++++++--------
 1 files changed, 169 insertions(+), 32 deletions(-)

diff --git a/src/mini_blas.c b/src/mini_blas.c
index 3af36e5..bac3e22 100644
--- a/src/mini_blas.c
+++ b/src/mini_blas.c
@@ -1,8 +1,12 @@
 
 #include <stdlib.h>
+#include <stdio.h>
 #include <math.h>
+#include <time.h>
+#include <string.h>
+#include "mini_blas.h"
 
-void pm(int M, int N, double *A)
+void pm(int M, int N, float *A)
 {
     int i,j;
     for(i =0 ; i < M; ++i){
@@ -14,37 +18,16 @@
     printf("\n");
 }
 
-void gemm(int TA, int TB, int M, int N, int K, double ALPHA, 
-                    double *A, int lda, 
-                    double *B, int ldb,
-                    double BETA,
-                    double *C, int ldc)
+void gemm(int TA, int TB, int M, int N, int K, float ALPHA, 
+        float *A, int lda, 
+        float *B, int ldb,
+        float BETA,
+        float *C, int ldc)
 {
-    // Assume TA = 0, beta = 1 LULZ
-    int i,j,k;
-    if(TB && !TA){
-        for(i = 0; i < M; ++i){
-            for(j = 0; j < N; ++j){
-                register double sum = 0;
-                for(k = 0; k < K; ++k){
-                    sum += ALPHA*A[i*lda+k]*B[k+j*ldb];
-                }
-                C[i*ldc+j] += sum;
-            }
-        }
-    }else{
-        for(i = 0; i < M; ++i){
-            for(k = 0; k < K; ++k){
-                register double A_PART = ALPHA*A[i*lda+k];
-                for(j = 0; j < N; ++j){
-                    C[i*ldc+j] += A_PART*B[k*ldb+j];
-                }
-            }
-        }
-    }
+    gpu_gemm( TA,  TB,  M, N, K, ALPHA,A,lda, B, ldb,BETA,C,ldc);
 }
 
-void im2row(double *image, int h, int w, int c, int size, int stride, double *matrix)
+void im2row(float *image, int h, int w, int c, int size, int stride, float *matrix)
 {
     int i;
     int mc = c;
@@ -64,7 +47,7 @@
         matrix[i] = image[pc*h*w+ph*w+pw];
     }
 }
-void im2col(double *image, int h, int w, int c, int size, int stride, double *matrix)
+void im2col(float *image, int h, int w, int c, int size, int stride, float *matrix)
 {
     int b,p;
     int blocks = ((h-size)/stride+1)*((w-size)/stride+1);
@@ -84,9 +67,9 @@
 }
 
 //From Berkeley Vision's Caffe!
-void im2col_cpu(double* data_im, const int channels,
+void im2col_cpu(float* data_im, const int channels,
         const int height, const int width, const int ksize, const int stride,
-        double* data_col) 
+        float* data_col) 
 {
     int c,h,w;
     int height_col = (height - ksize) / stride + 1;
@@ -106,3 +89,157 @@
     }
 }
 
+void col2im_cpu(float* data_col, const int channels,
+        const int height, const int width, const int ksize, const int stride,
+        float* data_im) 
+{
+    int c,h,w;
+    int height_col = (height - ksize) / stride + 1;
+    int width_col = (width - ksize) / stride + 1;
+    int channels_col = channels * ksize * ksize;
+    for ( c = 0; c < channels_col; ++c) {
+        int w_offset = c % ksize;
+        int h_offset = (c / ksize) % ksize;
+        int c_im = c / ksize / ksize;
+        for ( h = 0; h < height_col; ++h) {
+            for ( w = 0; w < width_col; ++w) {
+                data_im[(c_im * height + h * stride + h_offset) * width
+                    + w * stride + w_offset]+= data_col[(c * height_col + h) * width_col + w];
+            }
+        }
+    }
+}
+
+float *random_matrix(int rows, int cols)
+{
+    int i;
+    float *m = calloc(rows*cols, sizeof(float));
+    for(i = 0; i < rows*cols; ++i){
+        m[i] = (float)rand()/RAND_MAX;
+    }
+    return m;
+}
+
+void time_random_matrix(int TA, int TB, int m, int k, int n)
+{
+    float *a;
+    if(!TA) a = random_matrix(m,k);
+    else a = random_matrix(k,m);
+    int lda = (!TA)?k:m;
+    float *b;
+    if(!TB) b = random_matrix(k,n);
+    else b = random_matrix(n,k);
+    int ldb = (!TB)?n:k;
+
+    float *c = random_matrix(m,n);
+    int i;
+    clock_t start = clock(), end;
+    for(i = 0; i<1000; ++i){
+        cpu_gemm(TA,TB,m,n,k,1,a,lda,b,ldb,1,c,n);
+    }
+    end = clock();
+    printf("Matrix Multiplication %dx%d * %dx%d, TA=%d, TB=%d: %lf ms\n",m,k,k,n, TA, TB, (float)(end-start)/CLOCKS_PER_SEC);
+    free(a);
+    free(b);
+    free(c);
+}
+
+void test_blas()
+{
+    time_random_matrix(0,0,100,100,100); 
+    time_random_matrix(1,0,100,100,100); 
+    time_random_matrix(0,1,100,100,100); 
+    time_random_matrix(1,1,100,100,100); 
+
+    time_random_matrix(0,0,1000,100,100); 
+    time_random_matrix(1,0,1000,100,100); 
+    time_random_matrix(0,1,1000,100,100); 
+    time_random_matrix(1,1,1000,100,100); 
+
+
+}
+
+void time_gpu_random_matrix(int TA, int TB, int m, int k, int n)
+{
+    float *a;
+    if(!TA) a = random_matrix(m,k);
+    else a = random_matrix(k,m);
+    int lda = (!TA)?k:m;
+    float *b;
+    if(!TB) b = random_matrix(k,n);
+    else b = random_matrix(n,k);
+    int ldb = (!TB)?n:k;
+
+    float *c = random_matrix(m,n);
+    int i;
+    clock_t start = clock(), end;
+    for(i = 0; i<1000; ++i){
+        gpu_gemm(TA,TB,m,n,k,1,a,lda,b,ldb,1,c,n);
+    }
+    end = clock();
+    printf("Matrix Multiplication %dx%d * %dx%d, TA=%d, TB=%d: %lf ms\n",m,k,k,n, TA, TB, (float)(end-start)/CLOCKS_PER_SEC);
+    free(a);
+    free(b);
+    free(c);
+}
+
+void test_gpu_accuracy(int TA, int TB, int m, int k, int n)
+{
+    srand(0);
+    float *a;
+    if(!TA) a = random_matrix(m,k);
+    else a = random_matrix(k,m);
+    int lda = (!TA)?k:m;
+    float *b;
+    if(!TB) b = random_matrix(k,n);
+    else b = random_matrix(n,k);
+    int ldb = (!TB)?n:k;
+
+    float *c = random_matrix(m,n);
+    float *c_gpu = random_matrix(m,n);
+    memset(c, 0, m*n*sizeof(float));
+    memset(c_gpu, 0, m*n*sizeof(float));
+    int i;
+        //pm(m,k,b);
+        gpu_gemm(TA,TB,m,n,k,1,a,lda,b,ldb,1,c_gpu,n);
+        //pm(m, n, c_gpu);
+        cpu_gemm(TA,TB,m,n,k,1,a,lda,b,ldb,1,c,n);
+        //pm(m, n, c);
+    double sse = 0;
+    for(i = 0; i < m*n; ++i) {
+        //printf("%f %f\n", c[i], c_gpu[i]);
+        sse += pow(c[i]-c_gpu[i], 2);
+    }
+    printf("Matrix Multiplication %dx%d * %dx%d, TA=%d, TB=%d: %g MSE\n",m,k,k,n, TA, TB, sse/(m*n));
+    free(a);
+    free(b);
+    free(c);
+}
+
+void test_gpu_blas()
+{
+    test_gpu_accuracy(0,0,17,10,10); 
+    test_gpu_accuracy(1,0,17,10,10); 
+    test_gpu_accuracy(0,1,17,10,10); 
+    test_gpu_accuracy(1,1,17,10,10); 
+
+    test_gpu_accuracy(0,0,1000,10,100); 
+    test_gpu_accuracy(1,0,1000,10,100); 
+    test_gpu_accuracy(0,1,1000,10,100); 
+    test_gpu_accuracy(1,1,1000,10,100); 
+
+    time_gpu_random_matrix(0,0,1000,1000,100); 
+    time_random_matrix(0,0,1000,1000,100); 
+
+    time_gpu_random_matrix(0,1,1000,1000,100); 
+    time_random_matrix(0,1,1000,1000,100); 
+
+    time_gpu_random_matrix(1,0,1000,1000,100); 
+    time_random_matrix(1,0,1000,1000,100); 
+
+    time_gpu_random_matrix(1,1,1000,1000,100); 
+    time_random_matrix(1,1,1000,1000,100); 
+
+}
+
+

--
Gitblit v1.10.0