cleaned up driver
[RBC.git] / driver.cu
1 /* This file is part of the Random Ball Cover (RBC) library.
2 * (C) Copyright 2010, Lawrence Cayton [lcayton@tuebingen.mpg.de]
3 */
4
5 #include<stdio.h>
6 #include<stdlib.h>
7 #include<cuda.h>
8 #include<sys/time.h>
9 #include<math.h>
10 #include "defs.h"
11 #include "utils.h"
12 #include "utilsGPU.h"
13 #include "rbc.h"
14 #include "brute.h"
15 #include "sKernel.h"
16
17 void parseInput(int,char**);
18 void readData(char*,matrix);
19 void readDataText(char*,matrix);
20 void evalNNerror(matrix, matrix, unint*);
21 void evalKNNerror(matrix,matrix,intMatrix);
22
23 char *dataFileX, *dataFileQ, *outFile;
24 char runBrute=0, runEval=0;
25 unint n=0, m=0, d=0, numReps=0;
26 unint deviceNum=0;
27 int main(int argc, char**argv){
28 matrix x, q;
29 intMatrix nnsRBC;
30 matrix distsRBC;
31 struct timeval tvB,tvE;
32 cudaError_t cE;
33 rbcStruct rbcS;
34
35 printf("*****************\n");
36 printf("RANDOM BALL COVER\n");
37 printf("*****************\n");
38
39 parseInput(argc,argv);
40
41 printf("Using GPU #%d\n",deviceNum);
42 if(cudaSetDevice(deviceNum) != cudaSuccess){
43 printf("Unable to select device %d.. exiting. \n",deviceNum);
44 exit(1);
45 }
46
47 size_t memFree, memTot;
48 cudaMemGetInfo(&memFree, &memTot);
49 printf("GPU memory free = %lu/%lu (MB) \n",(unsigned long)memFree/(1024*1024),(unsigned long)memTot/(1024*1024));
50
51 initMat( &x, n, d );
52 initMat( &q, m, d );
53 x.mat = (real*)calloc( sizeOfMat(x), sizeof(*(x.mat)) );
54 q.mat = (real*)calloc( sizeOfMat(q), sizeof(*(q.mat)) );
55
56 //Load data
57 readData( dataFileX, x );
58 readData( dataFileQ, q );
59
60 //Allocate space for NNs and dists
61 initIntMat( &nnsRBC, m, K );
62 initMat( &distsRBC, m, K );
63 nnsRBC.mat = (unint*)calloc( sizeOfIntMat(nnsRBC), sizeof(*nnsRBC.mat) );
64 distsRBC.mat = (real*)calloc( sizeOfMat(distsRBC), sizeof(*distsRBC.mat) );
65
66 printf("\nrunning rbc..\n");
67 //Build the RBC
68 gettimeofday(&tvB,NULL);
69 buildRBC(x, &rbcS, numReps, numReps);
70 gettimeofday(&tvE,NULL);
71 printf("\t.. build time for rbc = %6.4f \n",timeDiff(tvB,tvE));
72
73 //This finds the 32-NNs; if you are only interested in the 1-NN, use queryRBC(..) instead
74 gettimeofday(&tvB,NULL);
75 kqueryRBC(q, rbcS, nnsRBC, distsRBC);
76 gettimeofday(&tvE,NULL);
77 printf("\t.. query time for krbc = %6.4f \n",timeDiff(tvB,tvE));
78
79 if( runBrute ){
80 intMatrix nnsBrute;
81 matrix distsBrute;
82 initIntMat( &nnsBrute, m, K );
83 nnsBrute.mat = (unint*)calloc( sizeOfIntMat(nnsBrute), sizeof(*nnsBrute.mat) );
84 initMat( &distsBrute, m, K );
85 distsBrute.mat = (real*)calloc( sizeOfMat(distsBrute), sizeof(*distsBrute.mat) );
86
87 printf("running k-brute force..\n");
88 gettimeofday(&tvB,NULL);
89 bruteK(x,q,nnsBrute,distsBrute);
90 gettimeofday(&tvE,NULL);
91 printf("\t.. time elapsed = %6.4f \n",timeDiff(tvB,tvE));
92
93 free(nnsBrute.mat);
94 free(distsBrute.mat);
95 }
96
97 cE = cudaGetLastError();
98 if( cE != cudaSuccess ){
99 printf("Execution failed; error type: %s \n", cudaGetErrorString(cE) );
100 }
101
102 if( runEval )
103 evalKNNerror(x,q,nnsRBC);
104
105 destroyRBC(&rbcS);
106 cudaThreadExit();
107 free(nnsRBC.mat);
108 free(distsRBC.mat);
109 free(x.mat);
110 free(q.mat);
111 }
112
113
114 void parseInput(int argc, char **argv){
115 int i=1;
116 if(argc <= 1){
117 printf("\nusage: \n testRBC -x datafileX -q datafileQ -n numPts (DB) -m numQueries -d dim -r numReps [-o outFile] [-g GPU num] [-b] [-e]\n\n");
118 printf("\tdatafileX = binary file containing the database\n");
119 printf("\tdatafileQ = binary file containing the queries\n");
120 printf("\tnumPts = size of database\n");
121 printf("\tnumQueries = number of queries\n");
122 printf("\tdim = dimensionailty\n");
123 printf("\tnumReps = number of representatives\n");
124 printf("\toutFile = output file (optional); stored in text format\n");
125 printf("\tGPU num = ID # of the GPU to use (optional) for multi-GPU machines\n");
126 printf("\n\tuse -b to run brute force in addition the RBC\n");
127 printf("\tuse -e option to run evaluation routine\n");
128 printf("\n\n");
129 exit(0);
130 }
131
132 while(i<argc){
133 if(!strcmp(argv[i], "-x"))
134 dataFileX = argv[++i];
135 else if(!strcmp(argv[i], "-q"))
136 dataFileQ = argv[++i];
137 else if(!strcmp(argv[i], "-n"))
138 n = atoi(argv[++i]);
139 else if(!strcmp(argv[i], "-m"))
140 m = atoi(argv[++i]);
141 else if(!strcmp(argv[i], "-d"))
142 d = atoi(argv[++i]);
143 else if(!strcmp(argv[i], "-r"))
144 numReps = atoi(argv[++i]);
145 else if(!strcmp(argv[i], "-o"))
146 outFile = argv[++i];
147 else if(!strcmp(argv[i], "-g"))
148 deviceNum = atoi(argv[++i]);
149 else{
150 fprintf(stderr,"%s : unrecognized option.. exiting\n",argv[i]);
151 exit(1);
152 }
153 i++;
154 }
155
156 if( !n || !m || !d || !numReps || !dataFileX || !dataFileQ ){
157 fprintf(stderr,"more arguments needed.. exiting\n");
158 exit(1);
159 }
160
161 if(numReps>n){
162 fprintf(stderr,"can't have more representatives than points.. exiting\n");
163 exit(1);
164 }
165 }
166
167
168 void readData(char *dataFile, matrix x){
169 unint i;
170 FILE *fp;
171 unint numRead;
172
173 fp = fopen(dataFile,"r");
174 if(fp==NULL){
175 fprintf(stderr,"error opening file.. exiting\n");
176 exit(1);
177 }
178
179 for( i=0; i<x.r; i++ ){ //can't load everything in one fread
180 //because matrix is padded.
181 numRead = fread( &x.mat[IDX( i, 0, x.ld )], sizeof(real), x.c, fp );
182 if(numRead != x.c){
183 fprintf(stderr,"error reading file.. exiting \n");
184 exit(1);
185 }
186 }
187 fclose(fp);
188 }
189
190
191 void readDataText(char *dataFile, matrix x){
192 FILE *fp;
193 real t;
194 int i,j;
195
196 fp = fopen(dataFile,"r");
197 if(fp==NULL){
198 fprintf(stderr,"error opening file.. exiting\n");
199 exit(1);
200 }
201
202 for(i=0; i<x.r; i++){
203 for(j=0; j<x.c; j++){
204 if(fscanf(fp,"%f ", &t)==EOF){
205 fprintf(stderr,"error reading file.. exiting \n");
206 exit(1);
207 }
208 x.mat[IDX( i, j, x.ld )]=(real)t;
209 }
210 }
211 fclose(fp);
212 }
213
214
215 //find the error rate of a set of NNs, then print it.
216 void evalNNerror(matrix x, matrix q, unint *NNs){
217 struct timeval tvB, tvE;
218 unint i;
219
220 printf("\nComputing error rates (this might take a while)\n");
221 real *ranges = (real*)calloc(q.pr,sizeof(*ranges));
222 for(i=0;i<q.r;i++){
223 if(NNs[i]>n) printf("error");
224 ranges[i] = distVec(q,x,i,NNs[i]) - 10e-6;
225 }
226
227 unint *cnts = (unint*)calloc(q.pr,sizeof(*cnts));
228 gettimeofday(&tvB,NULL);
229 bruteRangeCount(x,q,ranges,cnts);
230 gettimeofday(&tvE,NULL);
231
232 long int nc=0;
233 for(i=0;i<m;i++){
234 nc += cnts[i];
235 }
236 double mean = ((double)nc)/((double)m);
237 double var = 0.0;
238 for(i=0;i<m;i++) {
239 var += (((double)cnts[i])-mean)*(((double)cnts[i])-mean)/((double)m);
240 }
241 printf("\tavg rank = %6.4f; std dev = %6.4f \n\n", mean, sqrt(var));
242 printf("(range count took %6.4f) \n", timeDiff(tvB, tvE));
243
244 if(outFile){
245 FILE* fp = fopen(outFile, "a");
246 fprintf( fp, "%d %6.5f %6.5f \n", numReps, mean, sqrt(var) );
247 fclose(fp);
248 }
249
250 free(ranges);
251 free(cnts);
252 }
253
254
255 //evals the error rate of k-nns
256 void evalKNNerror(matrix x, matrix q, intMatrix NNs){
257 struct timeval tvB, tvE;
258 unint i,j,k;
259
260 unint m = q.r;
261 printf("\nComputing error rates (this might take a while)\n");
262
263 unint *ol = (unint*)calloc( q.r, sizeof(*ol) );
264
265 intMatrix NNsB;
266 NNsB.r=q.r; NNsB.pr=q.pr; NNsB.c=NNsB.pc=32; NNsB.ld=NNsB.pc;
267 NNsB.mat = (unint*)calloc( NNsB.pr*NNsB.pc, sizeof(*NNsB.mat) );
268 matrix distsBrute;
269 distsBrute.r=q.r; distsBrute.pr=q.pr; distsBrute.c=distsBrute.pc=K; distsBrute.ld=distsBrute.pc;
270 distsBrute.mat = (real*)calloc( distsBrute.pr*distsBrute.pc, sizeof(*distsBrute.mat) );
271
272 gettimeofday(&tvB,NULL);
273 bruteK(x,q,NNsB,distsBrute);
274 gettimeofday(&tvE,NULL);
275
276 //calc overlap
277 for(i=0; i<m; i++){
278 for(j=0; j<K; j++){
279 for(k=0; k<K; k++){
280 ol[i] += ( NNs.mat[IDX(i, j, NNs.ld)] == NNsB.mat[IDX(i, k, NNsB.ld)] );
281 }
282 }
283 }
284
285 long int nc=0;
286 for(i=0;i<m;i++){
287 nc += ol[i];
288 }
289
290 double mean = ((double)nc)/((double)m);
291 double var = 0.0;
292 for(i=0;i<m;i++) {
293 var += (((double)ol[i])-mean)*(((double)ol[i])-mean)/((double)m);
294 }
295 printf("\tavg overlap = %6.4f/%d; std dev = %6.4f \n", mean, K, sqrt(var));
296
297 FILE* fp;
298 if(outFile){
299 fp = fopen(outFile, "a");
300 fprintf( fp, "%d %6.5f %6.5f ", numReps, mean, sqrt(var) );
301 }
302
303 real *ranges = (real*)calloc(q.pr,sizeof(*ranges));
304 for(i=0;i<q.r;i++){
305 ranges[i] = distVec(q,x,i,NNs.mat[IDX(i, K-1, NNs.ld)]);
306 }
307
308 unint *cnts = (unint*)calloc(q.pr,sizeof(*cnts));
309 bruteRangeCount(x,q,ranges,cnts);
310
311 nc=0;
312 for(i=0;i<m;i++){
313 nc += cnts[i];
314 }
315 mean = ((double)nc)/((double)m);
316 var = 0.0;
317 for(i=0;i<m;i++) {
318 var += (((double)cnts[i])-mean)*(((double)cnts[i])-mean)/((double)m);
319 }
320 printf("\tavg actual rank of 32nd NN returned by the RBC = %6.4f; std dev = %6.4f \n\n", mean, sqrt(var));
321 printf("(brute k-nn took %6.4f) \n", timeDiff(tvB, tvE));
322
323 if(outFile){
324 fprintf( fp, "%6.5f %6.5f \n", mean, sqrt(var) );
325 fclose(fp);
326 }
327
328 free(cnts);
329 free(ol);
330 free(NNsB.mat);
331 free(distsBrute.mat);
332 }