load balancer with taking inti account load of objects
[charm.git] / src / ck-ldb / tm_mapping.c
1 #include <unistd.h>
2 #include <stdio.h>
3 #include <stdlib.h>
4 #include <string.h>
5 #include <float.h>
6 #include "treetimings.h"
7 #include "treematch.h"
8 #include <ctype.h>
9 #include <math.h>
10 #include <assert.h> 
11 #include "treemapping.h"
12  
13 #define TEST_ERROR(n) {if(n!=0){fprintf(stderr,"Error %d Line %d\n",n,__LINE__);exit(-1);}}
14 #undef DEBUG
15
16 #define LINE_SIZE 1000000
17
18 typedef struct{
19   int  val;
20   long key;
21 }hash_t;
22
23
24 typedef struct{
25   double val;
26   int key1;
27   int key2;
28 }hash2_t;
29
30
31 int nb_nodes(tm_topology_t *topology){
32   return topology->nb_nodes[topology->nb_levels-1];
33 }
34
35
36 void free_topology(tm_topology_t *topology){
37    int i;
38   for(i=0;i<topology->nb_levels;i++)
39     free(topology->node_id[i]);
40   free(topology->node_id);
41   free(topology->nb_nodes);
42   free(topology->arity);
43   free(topology);
44 }
45
46 double print_sol(int N,int *Value,double **comm, double **arch){
47   double sol;
48   int i,j;
49   double a,c;
50
51   sol=0;
52   for (i=0;i<N;i++){
53     for (j=i+1;j<N;j++){
54       c=comm[i][j];
55       a=arch[Value[i]][Value[j]];
56       //printf("T_%d_%d %f/%f=%f\n",i,j,c,a,c/a);
57       sol+=c/a;
58     }
59   }
60   for (i = 0; i < N; i++) {
61     printf("%d", Value[i]);
62     if(i<N-1)
63       printf(",");
64       
65   }
66   printf(" : %g\n",sol);
67
68   return sol;
69 }
70
71
72 void print_1D_tab(int *tab,int N){
73   int i;
74   for (i = 0; i < N; i++) {
75     printf("%d", tab[i]);
76     if(i<N-1)
77       printf(",");
78   }
79   printf("\n");
80
81 }
82
83 int nb_lines(char *filename){
84   FILE *pf;
85   char  line[LINE_SIZE];
86   int N=0;
87
88   if(!(pf=fopen(filename,"r"))){
89       fprintf(stderr,"Cannot open %s\n",filename);
90       exit(-1);
91   }
92
93   while(fgets(line,LINE_SIZE,pf))
94     N++;
95
96   printf("N=%d\n",N);
97
98   fclose(pf);
99   return N;
100 }
101
102 void init_comm(char *filename,int N,double **comm){
103   int i,j;
104   FILE *pf;
105   char *ptr;
106   char  line[LINE_SIZE];
107
108   if(!(pf=fopen(filename,"r"))){
109     fprintf(stderr,"Cannot open %s\n",filename);
110     exit(-1);
111   }
112
113
114   j=-1;
115   i=0;
116   while(fgets(line,LINE_SIZE,pf)){
117
118   
119     char *l=line;
120     j=0;
121     //printf("%s|",line);
122     while((ptr=strsep(&l," \t"))){
123       if((ptr[0]!='\n')&&(!isspace(ptr[0]))&&(*ptr)){
124         comm[i][j]=atof(ptr);
125         //printf ("comm[%d][%d]=%f|%s|\n",i,j,comm[i][j],ptr);
126         j++;
127       }
128     }
129     if(j!=N){
130       fprintf(stderr,"Error at %d %d (%d!=%d)for %s\n",i,j,j,N,filename);
131       exit(-1);
132     }
133     i++;
134   }
135   if(i!=N){
136     fprintf(stderr,"Error at %d %d for %s\n",i,j,filename);
137     exit(-1);
138   }
139
140   /*  printf("%s:\n",filename);
141   for(i=0;i<N;i++){
142     for(j=0;j<N;j++){
143       printf("%6.1f ",comm[i][j]);
144     }
145     printf("\n");
146     } */
147
148 }    
149
150 int  build_comm(char *filename,double ***pcomm){
151   int i;
152   int N;
153   double **comm;
154   N=nb_lines(filename);
155   comm=(double**)malloc(N*sizeof(double*));
156   for(i=0;i<N;i++)
157     comm[i]=(double*)malloc(N*sizeof(double));
158   init_comm(filename,N,comm);
159   *pcomm=comm;
160   return N;
161 }
162
163 void map_Packed(tm_topology_t *topology,int N,int *Value){
164   int i,depth;
165
166   depth=topology->nb_levels-1;
167
168   for(i=0;i<N;i++){
169     //printf ("%d -> %d\n",objs[i]->os_index,i);
170     Value[i]=topology->node_id[depth][i];
171   }
172
173 }
174
175 void map_RR(int N,int *Value){
176   int i;
177
178   for(i=0;i<N;i++){
179     //printf ("%d -> %d\n",i,i);
180     Value[i]=i;
181   }
182
183 }
184
185
186
187 int hash_asc(const void* x1,const void* x2){ 
188
189   hash_t *e1,*e2;
190
191   e1=((hash_t*)x1);
192   e2=((hash_t*)x2);
193
194   
195   return e1->key<e2->key?-1:1;
196
197
198
199 int *generate_random_sol(tm_topology_t *topology,int N,int level,int seed){
200   hash_t *hash_tab;
201   int *sol,i;
202   int *nodes_id;
203
204   nodes_id=topology->node_id[level];
205
206
207   hash_tab=(hash_t*)malloc(sizeof(hash_t)*N);
208   sol=(int*)malloc(sizeof(int)*N);
209   
210   srandom(seed);
211   
212   for(i=0;i<N;i++){
213     hash_tab[i].val=nodes_id[i];
214     hash_tab[i].key=random();
215   }
216   
217   qsort(hash_tab,N,sizeof(hash_t),hash_asc);
218   for(i=0;i<N;i++){
219     sol[i]=hash_tab[i].val;
220   }
221   free(hash_tab);
222   return sol;
223 }
224
225
226 double eval_sol(int *sol,int N,double **comm, double **arch){
227   double res;
228   int i,j;
229   double a,c;
230
231   res=0;
232   for (i=0;i<N;i++){
233     for (j=i+1;j<N;j++){
234       c=comm[i][j];
235       a=arch[sol[i]][sol[j]];
236       res+=c/a;
237     }
238   }
239   return res;
240 }
241
242 void exchange(int *sol,int i,int j){
243   int tmp;
244   tmp=sol[i];
245   sol[i]=sol[j];
246   sol[j]=tmp;
247 }
248
249 double gain_exchange(int *sol,int l,int m,double eval1,int N,double **comm, double **arch){
250   double eval2;
251   if(l==m)
252     return 0;
253   exchange(sol,l,m);
254   eval2=eval_sol(sol,N,comm,arch);
255   exchange(sol,l,m);
256   return eval1-eval2;
257 }
258
259 void select_max(int *l,int *m,double **gain,int N,int *state){
260   int i,j;
261   double max;
262   max=-DBL_MAX;
263   
264   for(i=0;i<N;i++){
265     if(!state[i]){
266       for(j=0;j<N;j++){
267           if((i!=j)&&(!state[j])){
268             if(gain[i][j]>max){
269               *l=i;*m=j;
270               max=gain[i][j];
271             }
272           }
273       }
274     }
275   }
276   
277 }
278
279 void compute_gain(int *sol,int N,double **gain,double **comm, double **arch){
280   int i,j;
281   double eval1;
282   eval1=eval_sol(sol,N,comm,arch);
283   for(i=0;i<N;i++){
284     for(j=0;j<=i;j++){
285       gain[i][j]=gain[j][i]=gain_exchange(sol,i,j,eval1,N,comm,arch);
286     }
287   } 
288 }
289
290
291
292 /* Randomized Algorithm of
293 Hu Chen, Wenguang Chen, Jian Huang ,Bob Robert,and H.Kuhn. Mpipp: an automatic profile-guided
294 parallel process placement toolset for smp clusters and multiclusters. In
295 Gregory K. Egan and Yoichi Muraoka, editors, ICS, pages 353-360. ACM, 2006.
296  */
297
298 void map_MPIPP(tm_topology_t *topology,int nb_seed,int N,int *Value,double **comm, double **arch){
299   int *sol;
300   int *state;
301   double **gain;
302   int **history;
303   double *temp;
304   int i,j,t,l=0,m=0,loop=0,seed=0;
305   double max,sum,best_eval,eval;
306
307   
308   gain=(double**)malloc(sizeof(double*)*N);
309   for(i=0;i<N;i++){
310     gain[i]=(double*)malloc(sizeof(double)*N);  
311     if(!gain[i]){
312     }
313   }
314   history=(int**)malloc(sizeof(int*)*N);
315   for(i=0;i<N;i++)
316     history[i]=(int*)malloc(sizeof(int)*3);
317
318   state=(int*)malloc(sizeof(int)*N);
319   temp=(double*)malloc(sizeof(double)*N);
320
321   sol=generate_random_sol(topology,N,topology->nb_levels-1,seed++);
322   for(i=0;i<N;i++)
323     Value[i]=sol[i];
324   
325   best_eval=DBL_MAX;
326   while(seed<=nb_seed){
327     loop=0;
328     do{
329
330       for(i=0;i<N;i++){
331         state[i]=0;
332         //printf("%d ",sol[i]);
333       }
334       //printf("\n");
335       compute_gain(sol,N,gain,comm,arch);
336   
337       //display_tab(gain,N);
338       //exit(-1);
339       for(i=0;i<N/2;i++){
340         select_max(&l,&m,gain,N,state);
341         //printf("%d: %d <=> %d : %f\n",i,l,m,gain[l][m]);
342         state[l]=1;state[m]=1;
343         exchange(sol,l,m);
344         history[i][1]=l;history[i][2]=m;
345         temp[i]=gain[l][m];
346         compute_gain(sol,N,gain,comm,arch);
347       }
348
349       t=-1;
350       max=0;
351       sum=0;
352       for(i=0;i<N/2;i++){
353         sum+=temp[i];
354         if(sum>max){
355           max=sum;
356           t=i;
357         }
358       }
359       /*for(j=0;j<=t;j++)
360         printf("exchanging: %d with %d for gain: %f\n",history[j][1],history[j][2],temp[j]); */
361       for(j=t+1;j<N/2;j++){
362         exchange(sol,history[j][1],history[j][2]);
363         //printf("Undoing: %d with %d for gain: %f\n",history[j][1],history[j][2],temp[j]); 
364       }
365       //printf("max=%f\n",max);
366
367       /*for(i=0;i<N;i++){
368         printf("%d ",sol[i]);
369         }
370         printf("\n");*/
371       eval=eval_sol(sol,N,comm,arch);
372       if(eval<best_eval){
373         best_eval=eval;
374         for(i=0;i<N;i++)
375           Value[i]=sol[i];
376         //print_sol(N);
377       }
378     
379
380     }while(max>0);
381     
382     sol=generate_random_sol(topology,N,topology->nb_levels-1,seed++);
383
384   }
385
386 }
387   
388
389
390
391 void map_tree(tree_t* t1,tree_t *t2){
392   /*  double x1,x2;
393   if((!t1->left)&&(!t1->right)){
394     printf ("%d -> %d\n",t1->id,t2->id);
395     Value[t2->id]=t1->id;
396    return;
397   }
398   x1=t2->right->val/t1->right->val+t2->left->val/t1->left->val;
399   x2=t2->left->val/t1->right->val+t2->right->val/t1->left->val;
400   if(x1<x2){
401     map_tree(t1->left,t2->left);
402     map_tree(t1->right,t2->right);
403   }else{
404     map_tree(t1->right,t2->left);
405     map_tree(t1->left,t2->right);
406     }*/
407 }
408
409 void depth_first(tree_t *comm_tree, int *proc_list,int *i){
410   int j;
411   if(!comm_tree->child){
412     proc_list[(*i)++]=comm_tree->id;
413     return;
414   }
415
416   for(j=0;j<comm_tree->arity;j++){
417     depth_first(comm_tree->child[j],proc_list,i);
418   }
419 }
420
421 int nb_leaves(tree_t *comm_tree){
422   int n=0,j;
423
424   if(!comm_tree->child){
425     return 1;
426   }
427
428   for(j=0;j<comm_tree->arity;j++){
429     n+=nb_leaves(comm_tree->child[j]);
430   }
431   return n;
432 }
433
434
435
436
437 /*Map topology to cores: 
438  sigma_i is such that  process i is mapped on core sigma_i
439  k_i is such that core i exectutes process k_i
440
441  size of sigma is the number of process
442  size of k is the number of cores/nodes
443
444  We must have numbe of process<=number of cores
445
446  k_i =-1 if no process is mapped on core i
447 */
448 void map_topology(tm_topology_t *topology,tree_t *comm_tree,int nb_proc,int level,
449                   int *sigma, int *k){
450   int *nodes_id;
451   int N;
452   int *proc_list,i,l;
453   int M;
454   int block_size;
455
456  
457   M=nb_leaves(comm_tree);
458   printf("nb_leaves=%d\n",M);
459
460
461
462   nodes_id=topology->node_id[level];
463   N=topology->nb_nodes[level];
464   //printf("level=%d, nodes_id=%p, N=%d\n",level,nodes_id,N);
465
466
467   //printf("N=%d,nb_proc=%d\n",N,nb_proc);
468   /* The number of node at level "level" in the tree should be equal to the number of processors*/
469   assert(N==nb_proc);
470
471
472   proc_list=(int*)malloc(sizeof(int)*M);
473   i=0;
474   depth_first(comm_tree,proc_list,&i);
475
476   l=0;
477   for(i=0;i<M;i++){
478     //printf ("%d\n",proc_list[i]);
479   }
480
481
482   block_size=M/N;
483
484
485   if(k){/*if we need the k vector*/
486     printf("M=%d, N=%d, BS=%d\n",M,N,block_size);
487     for(i=0;i<nb_nodes(topology);i++){
488       k[i]=-1;
489     }
490     for(i=0;i<M;i++){
491       if(proc_list[i]!=-1){
492 #ifdef DEBUG
493         printf ("%d->%d\n",proc_list[i],nodes_id[i/block_size]);
494 #endif
495         sigma[proc_list[i]]=nodes_id[i/block_size];
496         k[nodes_id[i/block_size]]=proc_list[i];
497       }
498     }
499   }else{
500     printf("M=%d, N=%d, BS=%d\n",M,N,block_size);
501     for(i=0;i<M;i++){
502       if(proc_list[i]!=-1){
503 #ifdef DEBUG
504         printf ("%d->%d\n",proc_list[i],nodes_id[i/block_size]);
505 #endif
506         sigma[proc_list[i]]=nodes_id[i/block_size];
507       }
508     }
509
510   }
511   free(proc_list);
512
513 }
514
515 void map_topology_simple(tm_topology_t *topology,tree_t *comm_tree, int *sigma,int *k){
516   map_topology(topology,comm_tree,topology->nb_nodes[topology->nb_levels-1],topology->nb_levels-1,sigma,k);
517 }
518
519 int int_cmp(const void* x1,const void* x2){ 
520
521   int *e1,*e2;
522
523   e1=((int *)x1);
524   e2=((int *)x2);
525
526   
527   return (*e1)>(*e2)?-1:1;
528
529
530
531 int decompose(int n,int optimize,int *tab){
532   int primes[6]={2,3,5,7,11,0},i=0;
533   int flag=2;
534   int j=1;
535
536
537   while(primes[i]&&(n!=1)){
538     printf("[%d] before=%d\n",primes[i],n);
539     if(flag&&optimize&&(n%primes[i]!=0)){
540       n+=primes[i]-n%primes[i];
541       flag--;
542       i=0;
543       continue;
544     }
545     printf("after=%d\n",n);
546     if(n%primes[i]==0){
547       tab[j++]=primes[i];
548       n/=primes[i];
549     }else{
550       i++;
551       flag=1;
552     }
553   }
554   if(n!=1){
555     tab[j++]=n;
556   }
557
558   qsort(tab+1,j-1,sizeof(int),int_cmp);
559
560   for(i=0;i<j;i++)
561     printf("%d:",tab[i]);
562   printf("\n");
563
564   tab[j]=0;
565   
566   return j+1;
567 }
568
569
570 tree_t *build_synthetic_topology(int *synt_tab,int id,int depth,int nb_levels){
571   tree_t *res,**child;
572   int arity=synt_tab[0];
573   int val,i; 
574
575   res=(tree_t*)malloc(sizeof(tree_t));
576   val=0;
577   if(depth>=nb_levels)
578     child=NULL;
579   else{
580     child=(tree_t**)malloc(sizeof(tree_t*)*arity);
581     for(i=0;i<arity;i++){
582       child[i]=build_synthetic_topology(synt_tab+1,i,depth+1,nb_levels);
583     child[i]->parent=res;
584     val+=child[i]->val;
585     }
586   }
587   set_node(res,child,arity,NULL,id,val+speed(depth),child[0]);
588   return res;
589 }
590
591
592 void   build_synthetic_proc_id(tm_topology_t *topology){
593   int n=1,i,j;
594   topology->node_id=(int**)malloc(sizeof(int*)*topology->nb_levels);
595   topology->nb_nodes=(int*)malloc(sizeof(int)*topology->nb_levels);
596
597
598   for(i=0;i<topology->nb_levels;i++){
599     topology->nb_nodes[i]=n;
600     topology->node_id[i]=(int*)malloc(sizeof(int)*n);
601     for(j=0;j<n;j++)
602       topology->node_id[i][j]=j;
603     n*=topology->arity[i]; 
604   }
605
606 }
607
608 void TreeMatchMapping(int nb_obj, int nb_proc,double **comm_mat,  int *sol){
609   tree_t *comm_tree;
610   tm_topology_t *topology;
611
612   int i;
613   for(i=0;i<nb_obj;i++)
614     sol[i]=i;
615
616
617   //  return;
618
619   topology=(tm_topology_t*)malloc(sizeof(tm_topology_t));
620   topology->arity=(int*)malloc(sizeof(int)*MAX_LEVELS);
621   topology->arity[0]=nb_proc;
622   topology->nb_levels=decompose((int)ceil((1.0*nb_obj)/nb_proc),1,topology->arity);
623   printf("Topology nb levels=%d\n",topology->nb_levels);
624   build_synthetic_proc_id(topology);
625
626
627   //exit(-1);
628   //topology_to_arch(topology);
629
630   //display_tab(arch,hwloc_get_nbobjs_by_type(topology, HWLOC_OBJ_PROC));
631   //display_tab(arch,96);
632   //exit(-1);
633   //int nb_core=topo_nb_proc(topology,1000);
634
635   //display_tab(comm_mat,N);
636
637   comm_tree=build_tree_from_topology(topology,comm_mat,nb_obj);
638   map_topology(topology,comm_tree,nb_proc,1,sol,NULL);
639
640
641   free_topology(topology);
642   free_tree(comm_tree);
643   printf("-------------- Mapping done!\n");
644 }
645
646 void display_other_heuristics(tm_topology_t *topology,int N,double **comm,double **arch){
647   CLOCK_T time1,time0;
648   double duration; 
649   int *sol;
650  
651   sol=(int*)malloc(sizeof(int)*N);
652
653   map_Packed(topology,N,sol);
654   printf("Packed: "); 
655   print_sol(N,sol,comm,arch);
656
657
658
659   map_RR(N,sol);
660   printf("RR: "); 
661   print_sol(N,sol,comm,arch);
662
663   CLOCK(time0);
664   map_MPIPP(topology,1,N,sol,comm,arch);
665   CLOCK(time1);
666   duration=CLOCK_DIFF(time1,time0);
667   printf("MPIPP-1-D:%f\n",duration);
668   printf("MPIPP-1: ");
669   print_sol(N,sol,comm,arch);
670
671   CLOCK(time0);
672   map_MPIPP(topology,5,N,sol,comm,arch);
673   CLOCK(time1);
674   duration=CLOCK_DIFF(time1,time0);
675   printf("MPIPP-5-D:%f\n",duration);
676   printf("MPIPP-5: ");
677   print_sol(N,sol,comm,arch);
678
679   free(sol);
680 }
681
682