#define NORMALIZE_T 1#define REINELT_START 0 \#define MAX_ITERS 9999 \/*6:*/#line 249 "./ascend.w"#include <config.h>#include "lkconfig.h"/*7:*/#line 263 "./ascend.w"#include <stdio.h>#include <stdlib.h>#include <stddef.h>/*:7*/#line 252 "./ascend.w"/*8:*/#line 269 "./ascend.w"#include "error.h"#include "length.h"#include "read.h" #include "lk.h"#include "memory.h"#include "decluster.h"#include "dict.h"#include "pq.h"#include "nn.h"/*:8*//*10:*/#line 292 "./ascend.w"#include "ascend.h"/*:10*/#line 253 "./ascend.w"/*15:*/#line 349 "./ascend.w"static int n= 0;static int total_iter= 0;static double*best_lambda= NULL,*read_lambda= NULL,*write_lambda= NULL;/*:15*//*18:*/#line 370 "./ascend.w"static int*degree_less_2= NULL;static int*degree_less_2_old= NULL;/*:18*//*30:*/#line 574 "./ascend.w"static decluster_edge_t*onetree= NULL;/*:30*//*46:*/#line 780 "./ascend.w"static int*work_from= NULL;static length_t*work_dist= NULL;/*:46*//*52:*/#line 860 "./ascend.w"static dict_t*master_edges= NULL;static pool_t*edge_pool= NULL;static pq_t*edges_pq= NULL;/*:52*//*64:*/#line 1010 "./ascend.w"int*parent,*count;/*:64*/#line 255 "./ascend.w"/*20:*/#line 397 "./ascend.w"static double normalizer(int n,int*degree_less_2);/*:20*//*22:*/#line 424 "./ascend.w"static void update_lambda(int n,double t,int*degree_less_2,int*degree_less_2_old,double*read_lambda,double*write_lambda);/*:22*//*42:*/#line 745 "./ascend.w"static length_t compute_onetree(decluster_edge_t*onetree,int approximate_mst);/*:42*//*56:*/#line 923 "./ascend.w"void add_edges_from(decluster_tree_t*mst);/*:56*//*60:*/#line 964 "./ascend.w"static void write_edge_weight_and_insert(void*env2,void**payload_p);/*:60*/#line 256 "./ascend.w"/*47:*/#line 788 "./ascend.w"static length_tcustom_cost(int i,int j){return cost(i,j)+read_lambda[i]+read_lambda[j];}/*:47*//*61:*/#line 968 "./ascend.w"static length_tkruskal(int n,decluster_tree_t*T,pq_t*edges_pq){length_t len= 0;int i,j,num_edges,root[2],u,v;for(i= 0;i<n;i++){parent[i]= -1;count[i]= 1;}for(num_edges= 0;num_edges<n-1;){decluster_edge_t*e= pq_delete_min(edges_pq);for(j= 0;j<2;j++){for(u= e->city[j];parent[u]>=0;u= parent[u]);root[j]= u;for(v= e->city[j];v!=u;v= parent[v])parent[v]= u;}if(root[0]!=root[1]){T->edge[num_edges++]= *e;len+= e->cost;if(count[root[0]]>count[root[1]]){parent[root[1]]= root[0];count[root[0]]+= count[root[1]];}else{parent[root[0]]= root[1];count[root[1]]+= count[root[0]];}}}return len;}/*:61*/#line 257 "./ascend.w"/*19:*/#line 384 "./ascend.w"static doublenormalizer(int n,int*degree_less_2){int i;double l;for(l= 0.0,i= 0;i<n;i++){l+= degree_less_2[i]*degree_less_2[i];}return l;}/*:19*//*21:*/#line 409 "./ascend.w"static voidupdate_lambda(int n,double t,int*degree_less_2,int*degree_less_2_old,double*read_lambda,double*write_lambda){int i;const double recent_share= 0.75;const double ta= t*recent_share,tb= t*(1.0-recent_share);for(i= 0;i<n;i++){write_lambda[i]= read_lambda[i]+ta*degree_less_2[i]+tb*degree_less_2_old[i];}}/*:21*//*41:*/#line 734 "./ascend.w"static length_tcompute_onetree(decluster_edge_t*onetree,int approximate_mst){length_t len;/*43:*/#line 751 "./ascend.w"{decluster_tree_t T;T.n= n-2;T.edge= onetree;if(approximate_mst){/*58:*/#line 944 "./ascend.w"if(edges_pq){pq_make_empty(edges_pq);}else{edges_pq= pq_create_size(decluster_edge_cmp,(int)dict_size(master_edges));}dict_update_all(master_edges,write_edge_weight_and_insert,edges_pq);/*:58*/#line 757 "./ascend.w"len= kruskal(n-1,&T,edges_pq);}else{len= decluster_mst_custom(&T,work_from,work_dist,custom_cost);add_edges_from(&T);}}/*:43*/#line 739 "./ascend.w"/*48:*/#line 801 "./ascend.w"{const int v= n-1;int i,short_to[2];length_t short_dist[2]= {INFINITY,INFINITY};for(i= 0;i<v;i++){const length_t di= custom_cost(v,i);if(di<short_dist[0]){short_to[1]= short_to[0];short_dist[1]= short_dist[0];short_to[0]= i;short_dist[0]= di;}else if(di<short_dist[1]){short_to[1]= i;short_dist[1]= di;}}onetree[n-2].city[0]= v;onetree[n-2].city[1]= short_to[0];onetree[n-2].cost= short_dist[0];onetree[n-1].city[0]= v;onetree[n-1].city[1]= short_to[1];onetree[n-1].cost= short_dist[1];len+= short_dist[0]+short_dist[1];}/*:48*/#line 740 "./ascend.w"return len;}/*:41*//*55:*/#line 906 "./ascend.w"voidadd_edges_from(decluster_tree_t*mst){int present= 0,i;decluster_edge_t*edge= pool_alloc(edge_pool);printf("add_edges_from: dict size before %d",dict_size(master_edges));for(i= 0;i<mst->n;i++){const int from= mst->edge[i].city[0];const int to= mst->edge[i].city[1];/*54:*/#line 887 "./ascend.w"if(from==n-1||to==n-1){present= 1;}else{if(from<to){edge->city[0]= from;edge->city[1]= to;}else{edge->city[0]= to;edge->city[1]= from;}{present= dict_insert(master_edges,edge);if(!present){printf("!");edge= pool_alloc(edge_pool);}}}/*:54*/#line 916 "./ascend.w"}if(present){pool_free(edge_pool,edge);}printf(" after %d\n",dict_size(master_edges));}/*:55*//*57:*/#line 929 "./ascend.w"static intcmp_edge(const void*aa,const void*bb){const decluster_edge_t*a= (const decluster_edge_t*)aa;const decluster_edge_t*b= (const decluster_edge_t*)bb;const int a0= a->city[0],a1= a->city[1];const int b0= b->city[0],b1= b->city[1];if(a0<b0||(a0==b0&&a1<b1))return-1;if(a0==b0&&a1==b1)return 0;return 1;}/*:57*//*59:*/#line 953 "./ascend.w"static voidwrite_edge_weight_and_insert(void*env2,void**payload_p){pq_t*pq= env2;decluster_edge_t*e= *(decluster_edge_t**)payload_p;e->cost= custom_cost(e->city[0],e->city[1]);pq_insert(pq,e);}/*:59*/#line 258 "./ascend.w"/*11:*/#line 307 "./ascend.w"voidascend_setup(int the_n){n= the_n;best_lambda= new_arr_of(double,n);read_lambda= new_arr_of(double,n);write_lambda= new_arr_of(double,n);/*16:*/#line 360 "./ascend.w"degree_less_2= new_arr_of(int,n);degree_less_2_old= new_arr_of(int,n);/*:16*//*28:*/#line 566 "./ascend.w"onetree= new_arr_of(decluster_edge_t,n);/*:28*//*44:*/#line 769 "./ascend.w"work_from= new_arr_of(int,n-1);work_dist= new_arr_of(length_t,n-1);/*:44*//*50:*/#line 848 "./ascend.w"master_edges= dict_create(cmp_edge,NULL);edge_pool= pool_create(sizeof(decluster_edge_t),n+50);edges_pq= NULL;/*:50*//*62:*/#line 1000 "./ascend.w"parent= new_arr_of(int,n);count= new_arr_of(int,n);/*:62*/#line 315 "./ascend.w"total_iter= 0;}/*:11*//*12:*/#line 321 "./ascend.w"voidascend_cleanup(void){free_mem(best_lambda);free_mem(read_lambda);free_mem(write_lambda);mem_deduct(sizeof(double)*3*n);/*17:*/#line 365 "./ascend.w"free_mem(degree_less_2);mem_deduct(sizeof(int)*n);free_mem(degree_less_2_old);mem_deduct(sizeof(int)*n);/*:17*//*29:*/#line 570 "./ascend.w"free_mem(onetree);mem_deduct(n*sizeof(decluster_edge_t));/*:29*//*45:*/#line 774 "./ascend.w"free_mem(work_from);free_mem(work_dist);mem_deduct((n-1)*(sizeof(int)+sizeof(length_t)));/*:45*//*51:*/#line 854 "./ascend.w"dict_destroy(master_edges,NULL);pool_destroy(edge_pool);if(edges_pq)pq_destroy(edges_pq);/*:51*//*63:*/#line 1005 "./ascend.w"free_mem(parent);mem_deduct(sizeof(int)*n);free_mem(count);mem_deduct(sizeof(int)*n);/*:63*/#line 329 "./ascend.w"n= 0;total_iter= 0;}/*:12*//*13:*/#line 335 "./ascend.w"double*constascend_best_lambda(void){return best_lambda;}/*:13*//*24:*/#line 467 "./ascend.w"length_tascend_alpha_beta(const int n,length_t upper_bound_len,double alpha,double beta){extern int verbose;int i,best_iter,iter= 0,new_is_best;double t,step_scale= alpha,norm,best_lower_bound= 0.0,onetree_len,best_actual_lower_bound= 0.0;double reinelt_step_scale= 1.0,vj_step_scale= 1.0;int approximate_mst= 0;errorif(LENGTH_TYPE_IS_INTEGRAL,"Held-Karp lower bound computations require length_t to be a ""floating point type.  Sorry, but you have to recompile.");errorif(n<3,"ascend: n=%d < 3\n",n);errorif(beta>=1.0,"ascend: beta=%f > 1\n",beta);errorif(beta<=0,"ascend: beta=%f <= 0\n",beta);for(i= 0;i<n;i++)degree_less_2_old[i]= read_lambda[i]= 0.0;best_lower_bound= 0.0;/*53:*/#line 869 "./ascend.w"{int from,j,size,*neighbour,present= 0;extern decluster_tree_t*mst;decluster_edge_t*edge= pool_alloc(edge_pool);for(from= 0;from<n;from++){neighbour= nn_list(from,&size);for(j= 0;j<size;j++){const int to= neighbour[j];/*54:*/#line 887 "./ascend.w"if(from==n-1||to==n-1){present= 1;}else{if(from<to){edge->city[0]= from;edge->city[1]= to;}else{edge->city[0]= to;edge->city[1]= from;}{present= dict_insert(master_edges,edge);if(!present){printf("!");edge= pool_alloc(edge_pool);}}}/*:54*/#line 877 "./ascend.w"}}if(present){pool_free(edge_pool,edge);}add_edges_from(mst);}/*:53*/#line 486 "./ascend.w"while(1){approximate_mst= (iter%100>0);/*25:*/#line 528 "./ascend.w"/*31:*/#line 580 "./ascend.w"{const double len= compute_onetree(onetree,approximate_mst);double lambda_2= 0.0;int i;for(lambda_2= 0.0,i= 0;i<n;i++)lambda_2+= 2*read_lambda[i];onetree_len= len-lambda_2;if(verbose>=80&&!approximate_mst){printf("compute_onetree: A %s 1-tree is len %f\n",(approximate_mst?"approximate":"full"),onetree_len);fflush(stdout);}else{if(verbose){putchar('.');fflush(stdout);}}}/*:31*/#line 529 "./ascend.w"/*49:*/#line 829 "./ascend.w"{int i;for(i= 0;i<n;i++)degree_less_2[i]= -2;for(i= 0;i<n-1;i++){degree_less_2[onetree[i].city[0]]++;degree_less_2[onetree[i].city[1]]++;}}/*:49*/#line 530 "./ascend.w"if(!approximate_mst&&onetree_len>best_actual_lower_bound)best_actual_lower_bound= onetree_len;new_is_best= (onetree_len>best_lower_bound);norm= normalizer(n,degree_less_2);/*37:*/#line 664 "./ascend.w"if(new_is_best){best_lower_bound= onetree_len;best_iter= iter;}/*:37*/#line 535 "./ascend.w"/*:25*/#line 489 "./ascend.w"/*38:*/#line 672 "./ascend.w"if(new_is_best){double err= (upper_bound_len-onetree_len)/upper_bound_len;if(err<0.0005){if(verbose>=75)printf("#  Ascend: stopping criteria met: %.2f%% away from upper\n",err*100);break;}}/*:38*/#line 491 "./ascend.w"/*32:*/#line 605 "./ascend.w"if(norm==0.0){/*67:*/#line 1034 "./ascend.w"#if 0if(verbose>=50)printf("#  Ascend: Page %% total_iter %d Found a tour of length %f\n",total_iter,(float)(onetree_len));#endif/*:67*/#line 607 "./ascend.w"for(onetree_len= 0.0,i= 0;i<n;i++)onetree_len+= onetree[i].cost;best_lower_bound= onetree_len;/*68:*/#line 1043 "./ascend.w"#if 0if(verbose>=100){char s[100];sprintf(s,"Tour length %f",(float)onetree_len);show_onetree(debug_ps,s,n,NULL,edge);}#endif/*:68*/#line 611 "./ascend.w"break;}/*:32*/#line 492 "./ascend.w"/*35:*/#line 645 "./ascend.w"if(iter>=MAX_ITERS){if(verbose>=75)printf("#  Ascend: Iterations exceeded %d\n",MAX_ITERS);break;}/*:35*/#line 493 "./ascend.w"if(iter==0){reinelt_step_scale= 15*(upper_bound_len-onetree_len)/n;vj_step_scale= alpha;}else{reinelt_step_scale*= beta;vj_step_scale*= beta;}if(0&&(iter%500)==0){beta= 1-(1-beta)/2;}if(1||iter>100){step_scale= reinelt_step_scale;t= step_scale;}else{step_scale= vj_step_scale;t= step_scale*(upper_bound_len-onetree_len)/norm;}/*65:*/#line 1014 "./ascend.w"if(verbose>=100){printf("%d %f # L(l) t = %f beta=%f\n",total_iter,onetree_len,t,beta);fflush(stdout);}/*:65*/#line 510 "./ascend.w"/*36:*/#line 654 "./ascend.w"if(t<1e-3){if(verbose>=75)printf("#  Ascend: stopping criteria met: t < 0.001 \n");break;}/*:36*/#line 511 "./ascend.w"/*34:*/#line 620 "./ascend.w"#if 0if(is_random_dist_matrix){if(iter>1000&&iter-best_iter>5&&iter-best_iter<8){if(verbose>=75)printf("#  Ascend: stopping criteria met:  best is old\n");break;}}else{if(iter>1000&&iter-best_iter>5){if(verbose>=75)printf("#  Ascend: stopping criteria met:  best is old\n");break;}}#endif/*:34*/#line 512 "./ascend.w"update_lambda(n,t,degree_less_2,degree_less_2_old,read_lambda,write_lambda);/*66:*/#line 1022 "./ascend.w"#if 0if(total_iter==0&&verbose>=100)show_onetree(debug_ps,"First 1-tree, lambda==0 vector",n,NULL,edge);else if(verbose>=500)show_onetree(debug_ps,NULL,n,NULL,edge);if(verbose)fflush(stdout);#endif/*:66*/#line 515 "./ascend.w"/*39:*/#line 694 "./ascend.w"{double*b= best_lambda,*r= read_lambda,*w= write_lambda;if(new_is_best){best_lambda= b;read_lambda= w;write_lambda= r;}else{best_lambda= r;read_lambda= w;write_lambda= b;}}/*:39*/#line 516 "./ascend.w"/*23:*/#line 433 "./ascend.w"{int*t= degree_less_2;degree_less_2= degree_less_2_old;degree_less_2_old= t;}/*:23*/#line 517 "./ascend.w"total_iter++,iter++;}approximate_mst= 0;/*39:*/#line 694 "./ascend.w"{double*b= best_lambda,*r= read_lambda,*w= write_lambda;if(new_is_best){best_lambda= b;read_lambda= w;write_lambda= r;}else{best_lambda= r;read_lambda= w;write_lambda= b;}}/*:39*/#line 521 "./ascend.w"for(i= 0;i<n;i++)read_lambda[i]= best_lambda[i];/*25:*/#line 528 "./ascend.w"/*31:*/#line 580 "./ascend.w"{const double len= compute_onetree(onetree,approximate_mst);double lambda_2= 0.0;int i;for(lambda_2= 0.0,i= 0;i<n;i++)lambda_2+= 2*read_lambda[i];onetree_len= len-lambda_2;if(verbose>=80&&!approximate_mst){printf("compute_onetree: A %s 1-tree is len %f\n",(approximate_mst?"approximate":"full"),onetree_len);fflush(stdout);}else{if(verbose){putchar('.');fflush(stdout);}}}/*:31*/#line 529 "./ascend.w"/*49:*/#line 829 "./ascend.w"{int i;for(i= 0;i<n;i++)degree_less_2[i]= -2;for(i= 0;i<n-1;i++){degree_less_2[onetree[i].city[0]]++;degree_less_2[onetree[i].city[1]]++;}}/*:49*/#line 530 "./ascend.w"if(!approximate_mst&&onetree_len>best_actual_lower_bound)best_actual_lower_bound= onetree_len;new_is_best= (onetree_len>best_lower_bound);norm= normalizer(n,degree_less_2);/*37:*/#line 664 "./ascend.w"if(new_is_best){best_lower_bound= onetree_len;best_iter= iter;}/*:37*/#line 535 "./ascend.w"/*:25*/#line 523 "./ascend.w"return best_actual_lower_bound;}/*:24*//*26:*/#line 553 "./ascend.w"length_tascend(const int n,length_t upper_bound_len){return ascend_alpha_beta(n,upper_bound_len,1.5,0.999);}/*:26*/#line 259 "./ascend.w"const char*ascend_rcs_id= "$Id: ascend.w,v 1.37 1998/12/05 22:37:53 neto Exp neto $";/*:6*/