added utils, basic structure

cmake/blasJoin/blasJoin.cpp

+
+#include "blasJoin.h"
+
+double blasJoin(double *x, const size_t N, const size_t D, const double EPS, const unsigned int THREADS){
+    __builtin_assume_aligned(x, 64);
+    double elapsed=0.0;
+    CUtilTimer timer;
+    double *result;
+
+    if ( NUM_THREADS != 0 ){
+         omp_set_num_threads(THREADS);
+         mkl_set_num_threads(THREADS);
+    }
+
+    result = (double*) ddr_alloc(N * N);
+    // ah this shit does not fit into memory.
+    // obviously, N=200000 leads to a size of ~298GB. 
+    // TODO: make a boolean matrix as result, and do some shifting operation to
+    // determine wheater we have a positive or negavite number
+
+    // test BLAS:
+    timer.start();
+    cblas_dgemm(CblasRowMajor, CblasNoTrans, CblasTrans, N, N, D, 1.0, x, D, x, D, 0.0, result, N);
+    timer.stop();
+
+    ddr_free(result);
+    return timer.get_time();
+}

cmake/blasJoin/blasJoin.h

+
+#ifndef BLAS_JOIN_H
+#define BLAS_JOIN_H
+
+#include <stdio.h>
+#include <assert.h>
+#include <string.h>
+#include <omp.h>
+
+#include "mkl.h"
+
+#include "../util/allocation.h"
+#include "../util/timer.h"
+
+double blasJoin(double *x, const size_t N, const size_t D, const double EPS, const unsigned int THREADS);
+
+#endif

cmake/util/allocation.cpp

+
+#include "allocation.h"
+
+void * ddr_alloc(size_t bytes){
+    void * ptr = _mm_malloc(bytes, ALIGNMENT);
+
+    if ( ptr == NULL ){
+        fprintf(stderr, "Error in allocating memory with ddr_alloc!");
+        exit(1);
+    }
+
+    return ptr;
+}
+
+void ddr_free(void *ptrs){
+    _mm_free(ptrs);
+}

cmake/util/allocation.h

+
+#ifndef ALLOCATION_H
+#define ALLOCATION_H
+
+#include <stdio.h>
+#include <stdlib.h>
+
+#if defined(__INTEL_COMPILER)
+#include <malloc.h>
+#else
+#include <mm_malloc.h>
+#endif
+
+#define ALIGNMENT 64
+
+void * ddr_alloc(size_t bytes);
+void ddr_free(void * ptrs);
+
+#endif

cmake/util/arguments.cpp

+
+
+#include "arguments.h"
+
+
+void parsing_args(int argc, char* argv[], size_t *n, double *epsilon, size_t *d, size_t *threads, char *filename, bool isBinary){
+  char c;
+  FILE *file;
+
+  if ( argc < 4 ){
+    fprintf (stderr, "The parameters are obligatory.\n");
+    fprintf (stderr, "Usage: ./blasJoin ");
+
+    fprintf(stderr, "Obligatory parameters: \n");
+    fprintf(stderr, "n (number of objects in thousands (*1000))\ne (epsilon)\nd (dimensionality)\n");
+    fprintf(stderr, "Optional parameters: \n t number of threads\n\n");
+    fprintf(stderr, "f (filename) if there is no filename we use random generated data [0.0, 100.0)\n");
+    fprintf(stderr, "b use the -b argument without options to specify that it is a binary file.\n");
+    fprintf(stderr, "Example (with default values): ./blasMeans -n 200 -e 0.2 -d 20 -t 64\n");
+    exit(1);
+  }
+
+  while ((c = getopt(argc, argv, "bn:e:d:t:f:")) != -1) {
+	if ( optarg ){
+        switch(c){
+            case 'n':
+    		    *n = atol(optarg);
+    		    break;
+    		case 't':
+    		    *threads = atoi(optarg);
+    			break;
+            case 'd':
+                *d = atoi(optarg);
+                break;
+            case 'e':
+                *epsilon = atof(optarg);
+                break;
+            case 'f':
+                strcpy(filename, optarg);
+                break;
+            case 'b':
+                isBinary = true;
+                break;
+            case '?':
+              if (optopt == 'c')
+                fprintf (stderr, "Option -%c requires an argument.\n", optopt);
+              else if (isprint (optopt)){
+                  fprintf (stderr, "Unknown option `-%c'.\n", optopt);
+                  exit(1);                  
+              }
+              else
+                fprintf (stderr,
+                         "Unknown option character `\\x%x'.\n",
+                         optopt);
+              exit(1);
+    		default:
+    			break;
+        }
+	}
+  }
+
+}

cmake/util/arguments.h

+
+#ifndef MKM_ARGS_H
+#define MKM_ARGS_H
+
+#include <unistd.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <ctype.h>
+
+void parsing_args(int argc, char* argv[], size_t *n, double *epsilon, size_t *d, size_t *threads, char *filename, bool isBinary);
+
+#endif //KMEANS_ARGS_H

cmake/util/dataIo.cpp

+
+#include "dataIo.h"
+
+void random_init(double *array, const int N, const int D){
+    short unsigned seed[3];
+    int i;
+
+    seed[0]=1; seed[1]=1; seed[2]=2;
+
+    #pragma omp parallel for firstprivate(seed)
+    for ( i=0 ; i < N * D ; i++ ){
+        array[i] = erand48(seed) * 100.0;
+    }
+}
+
+void read_file(double *array, const int N, const int D, char filename[], bool isBinary){
+    FILE *fp;
+    size_t counts = 0;
+    size_t i=0,j=0;
+    char line[MAX_LINE_LENGTH];
+    char *token=NULL;
+    const char space[2] = " ";
+
+    fp = fopen(filename,"r");
+
+    if ( fp == NULL ){
+        fprintf(stderr, "File '%s' does not exists!", filename);
+        exit(1);
+    }
+
+    if ( isBinary ){
+        // read binary file, everything at once
+        counts = fread(array, sizeof(double) * N * D, 1, fp);
+
+        if ( counts == 0 ) {
+            fprintf(stderr, "Binary file '%s' could not be read. Wrong format.", filename);
+            exit(1);
+        }
+    }else{
+        // processing a text file
+        // format: there are D double values each line. Each value is separated by a space character.
+        // notice MAX_LINE_LENGTH = 2049
+        i = 0;
+        while ( fgets ( line, MAX_LINE_LENGTH, fp ) != NULL &&
+                i < N ) {
+            if ( line[0] != '%'){ // ignore '%' comment char
+                token = strtok(line, space);
+                j=0;
+
+                while ( token != NULL &&
+                        j < D ){
+                    array[i*D + j] = atof(token); // 0.0 if no valid conversion
+                    token = strtok(NULL, space);
+                    j++;
+                }
+                i++;
+            }
+        }
+    }
+
+    fclose(fp);
+}
+
+void save_binary_file(double *array, const int N, const int D, char filename[]){
+    FILE *fp=NULL;
+    size_t counts = 0;
+
+    fp = fopen(filename, "w");
+
+    if ( fp == NULL ){
+        fprintf(stderr, "Could not open file '%s'!", filename);
+        exit(1);
+    }
+
+    counts = fwrite(array,sizeof(double) * N * D, 1, fp);
+
+    if ( counts == 0 ){
+        fprintf(stderr, "Error in writing file '%s'. Abort.", filename);
+        exit(1);
+    }
+
+    fclose(fp);
+}
+
+void save_text_file(double *array, const int N, const int D, char filename[]){
+    FILE *fp=NULL;
+    size_t counts = 0;
+    size_t i=0, j=0;
+    char line[MAX_LINE_LENGTH];
+    char strDouble[50];
+
+    fp = fopen(filename, "w");
+
+    if ( fp == NULL ){
+        fprintf(stderr, "Could not open file '%s'!", filename);
+        exit(1);
+    }
+
+    for ( i=0 ; i < N ; i++ ){
+        strcpy(line, "");
+        for ( j=0 ; j < D ; j++ ){
+            strcpy(strDouble, "");
+            sprintf(strDouble, "%f ", array[i*D + j]);
+            strcat(line, strDouble);
+
+        }
+        fprintf(fp, "%s\n", line);
+
+    }
+
+    fclose(fp);
+}

cmake/util/dataIo.h

+
+#ifndef DATA_IO_H
+#define DATA_IO_H
+
+#include <stdlib.h>
+#include <string.h>
+#include <stdio.h>
+
+#include <omp.h>
+
+#define MAX_LINE_LENGTH 2049
+
+void random_init(double *array, const int N, const int D);
+void read_file(double *array, const int N, const int D, char filename[], bool isBinary);
+void save_binary_file(double *array, const int N, const int D, char filename[]);
+void save_text_file(double *array, const int N, const int D, char filename[]);
+
+#endif

cmake/util/timer.cpp

+//==============================================================
+//
+// SAMPLE SOURCE CODE - SUBJECT TO THE TERMS OF SAMPLE CODE LICENSE AGREEMENT,
+// http://software.intel.com/en-us/articles/intel-sample-source-code-license-agreement/
+//
+// Copyright 2013 Intel Corporation
+//
+// THIS FILE IS PROVIDED "AS IS" WITH NO WARRANTIES, EXPRESS OR IMPLIED, INCLUDING BUT
+// NOT LIMITED TO ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS FOR A PARTICULAR
+// PURPOSE, NON-INFRINGEMENT OF INTELLECTUAL PROPERTY RIGHTS.
+//
+// ===============================================================
+
+#include "timer.h"
+#ifdef _WIN32
+#include <Windows.h>
+#include <intrin.h>
+#else
+#include <sys/time.h>
+#endif
+
+// Description:
+// Registers the current clock tick value in m_start_clock_tick, current time value in m_start_time
+// Microsoft Windows* uses __rdtsc for clock ticks and QueryPerformanceFrequency/QueryPerformanceCounter for time
+// Linux*/OS X* uses the rdtsc instruction for clock ticks and get_timeofday for time
+void CUtilTimer::start() {
+#ifdef _WIN32
+	
+	// Clock ticks
+	//__rdtsc() is an intrinsic provided by Microsoft Visual Studio* in intrin.h header file
+	m_start_clock_tick = __rdtsc();
+
+	// Time
+	// QueryPerformanceFrequency works with QueryPerformanceCounter to return a human-readable time, provided in Windows.h
+	QueryPerformanceFrequency((LARGE_INTEGER *)&m_frequency);
+	unsigned __int64 now;
+	QueryPerformanceCounter((LARGE_INTEGER *)&now);
+	// Divide the raw counter by m_frequency for time in seconds
+	m_start_time = static_cast<double>(now) / m_frequency;
+#else
+
+	// Clock ticks
+	// On Linux and OS X, rdtsc instruction is used since we don't have intrinsic equivalent of __rdtsc()
+	unsigned lower, higher;
+	// rdtsc instruction returns a 64 bit clock tick 
+	// whose lower 32 bits is stored in EAX and higher 32 bits are stored in EDX register
+	__asm__ __volatile__("rdtsc":"=a"(lower), "=d"(higher));
+	// Constructing the 64 bit value from EAX and EDX
+	m_start_clock_tick = ((unsigned long long)lower)|(((unsigned long long)higher)<<32);
+
+	// Time
+	struct timeval start;
+	gettimeofday(&start, 0); //Returns the time of the day
+	//tv_sec records time in seconds and tv_usec records time in micro seconds
+	m_start_time = ((double) start.tv_sec + (double) start.tv_usec/1000000.0); 
+#endif
+}
+
+// Description:
+// Registers the current clock tick value in m_end_clock_tick, current time value in m_end_time
+// Windows uses __rdtsc for clock ticks and QueryPerformanceFrequency/QueryPerformanceCounter for time
+// Linux*/OS X* uses the rdtsc instruction for clock ticks and get_timeofday for time
+void CUtilTimer::stop() {
+#ifdef _WIN32
+
+	// Clock ticks
+	m_end_clock_tick = __rdtsc();
+
+	// Time
+	unsigned __int64 now;
+	QueryPerformanceCounter((LARGE_INTEGER *)&now);
+	m_end_time = static_cast<double>(now) / m_frequency;
+#else
+
+	// Clock ticks
+	unsigned lower, higher;
+	__asm__ __volatile__("rdtsc":"=a"(lower), "=d"(higher));
+	m_end_clock_tick = ((unsigned long long)lower)|(((unsigned long long)higher)<<32);
+
+	// Time
+	struct timeval start;
+	gettimeofday(&start, 0);
+	m_end_time = ((double) start.tv_sec + (double) start.tv_usec/1000000.0);
+#endif
+}
+
+// Description:
+// Returns the number of clock ticks taken between start and stop
+long long CUtilTimer::get_ticks() {
+	return (m_end_clock_tick - m_start_clock_tick);
+}
+
+// Description:
+// Returns the number of seconds taken between start and stop
+double CUtilTimer::get_time() {
+	return (m_end_time - m_start_time);
+}
+

cmake/util/timer.h

+//==============================================================
+//
+// SAMPLE SOURCE CODE - SUBJECT TO THE TERMS OF SAMPLE CODE LICENSE AGREEMENT,
+// http://software.intel.com/en-us/articles/intel-sample-source-code-license-agreement/
+//
+// Copyright 2013 Intel Corporation
+//
+// THIS FILE IS PROVIDED "AS IS" WITH NO WARRANTIES, EXPRESS OR IMPLIED, INCLUDING BUT
+// NOT LIMITED TO ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS FOR A PARTICULAR
+// PURPOSE, NON-INFRINGEMENT OF INTELLECTUAL PROPERTY RIGHTS.
+//
+// ===============================================================
+
+#ifndef TIMER_H
+#define TIMER_H
+
+class CUtilTimer {
+public:
+	CUtilTimer():
+		m_start_time(0.0),
+		m_end_time(0.0),
+		m_start_clock_tick(0),
+		m_end_clock_tick(0)
+	{};
+	// Registers the current clock tick and time value in m_start_clock_tick and m_start_time
+	void start();
+	// Registers the current clock tick and time value in m_end_clock_tick and m_end_time
+	void stop();
+	// Returns the number of seconds taken between start and stop
+	double get_time();
+	// Returns the number of clock ticks taken between start and stop
+	long long get_ticks();
+private:
+	// the start time and end time in seconds
+	double m_start_time, m_end_time; 
+	// the start clock tick and end clock tick
+	unsigned long long m_start_clock_tick, m_end_clock_tick;
+	// the frequency for QueryPerformance
+#ifdef _WIN32
+	unsigned __int64 m_frequency;
+#endif
+};
+
+#endif // TIMER_H