Want more features on Pastebin? Sign Up, it's FREE!
Guest

Untitled

By: a guest on Feb 26th, 2012  |  syntax: C++  |  size: 26.66 KB  |  views: 22  |  expires: Never
download  |  raw  |  embed  |  report abuse  |  print
Text below is selected. Please press Ctrl+C to copy to your clipboard. (⌘+C on Mac)
  1. #pragma warning(disable:4996)
  2. #include<stdio.h>
  3. #include<string.h>
  4. #include<math.h>
  5. #include<stdlib.h>
  6. #include<time.h>
  7. #include<algorithm>
  8. #include<string>
  9. #include<vector>
  10. #include<map>
  11. #include<set>
  12. #include<fstream>
  13.  
  14. #define inf 100000000
  15. #define MIN(a,b) ((a)<(b)?(a):(b))
  16. #define MAX(a,b) ((a)>(b)?(a):(b))
  17.  
  18. #define Tabu_prob 1                                 // rejects tabu listed substrings with "tabu_prob" probability
  19. #define Distinct_Substr 300005          // Maximum no. of distinct substrings
  20. #define NO_UPDATE 70                    // No update for NO_UPDATE iterations -> terminate
  21. #define MAX_RUN 4                              
  22. #define POPULATION_SIZE (n*2)           // No. of ants = some multiple of no. of target strings in the input
  23. #define T 305                                           // Maximum Cardinality of the set of target string
  24. #define M 305                                           // Maximum length of target string
  25. #define child 130                                       // no. of children of a node in TRIE
  26.  
  27. /* constants for a random number generator, for details see numerical recipes in C */
  28.  
  29. #define IA 16807
  30. #define IM 2147483647
  31. #define AM (1.0/IM)
  32. #define IQ 127773
  33. #define IR 2836
  34. #define MASK 123459876
  35. /* constants for a random number generator, for details see numerical recipes in C */
  36.  
  37. using namespace std;
  38.  
  39. // seed for ran01
  40. int seed = (int) time(NULL);
  41.  
  42. //is tabu on ? 1 : on, 0:Off
  43. int tabu_done;
  44.  
  45. char str_folder[100];
  46. double ran01( int *idum )
  47. /*    
  48.       FUNCTION:       generate a random number that is uniformly distributed in [0,1]
  49.       INPUT:          pointer to variable with the current seed
  50.       OUTPUT:         random number uniformly distributed in [0,1]
  51.       (SIDE)EFFECTS:  random number seed is modified (important, this has to be done!)
  52.       ORIGIN:         numerical recipes in C
  53. */
  54. {
  55.           long k;
  56.           double ans;
  57.  
  58.           k =(*idum)/IQ;
  59.           *idum = IA * (*idum - k * IQ) - IR * k;
  60.           if (*idum < 0 ) *idum += IM;
  61.           ans = AM * (*idum);
  62.           return ans;
  63. }
  64.  
  65. // v is target string idx, returns pair<> to access istarget
  66. pair<int,int> calc_target_idx(int v)
  67. {
  68.         return make_pair(v/30,v%30);
  69. }
  70.  
  71. struct trie
  72. {
  73. //      int freq_target;                // freq. of target strings
  74. //      int istarget[10];               // which target string contains this substring -> bitwise access
  75.         int val;                                // contains the unique idx of root to this position
  76. //      int freq;                               // total frequency of this substring
  77.         trie *next[child];
  78.         trie()
  79.         {
  80.                 int i;
  81.                 val=-1;
  82. //              memset(istarget,0,sizeof(istarget));
  83.                 for(i=0;i<child;i++)
  84.                         next[i]=NULL;
  85.         }
  86.         ~trie()
  87.         {
  88.                 int i;
  89.                 for(i=0;i<child;i++)
  90.                         delete next[i];
  91.         }
  92. };
  93.  
  94. //if q<=q0 then in ACS the best solution is selected
  95. const double q0 = 0.25;
  96. //prob of best
  97. const double p_best = 0.09;
  98. double t_max, t_min;                                    //Max, Min limit of pheromone value
  99. //this is not used -> will not work for MMAS
  100. const double minus_pheromone = 0.0;                                     //total amount of evaporation in a run
  101. //so it is made const. zero, so don't need to change anything in select_jump
  102. /********************************/
  103. double alpha, beta, evaporation_rate;  
  104. double allowed_time;                                            //may need to read from input
  105. int sub_freq[Distinct_Substr];                          //total freq. of substring w.r.t. it's id
  106. int sub_freq_target[Distinct_Substr];           //freq. of target strings who contain id-th substr
  107. int sub_is_target[Distinct_Substr][22];         //bitmap -> 1 -> which targets contain this id
  108. pair<int,int> sub_start[Distinct_Substr];
  109. pair<int,int> sub_end[Distinct_Substr];
  110. int sub_tabu[Distinct_Substr];                          //substring -> tabu or not
  111. int tabu_active;                                                        //marker of active tabu
  112. trie *root;
  113. int L;                                                                          //L of L-cover
  114. int id;                                                                         //gives id to a substring
  115. int n;                                                                          //number of target strings
  116. int m;                                                                          //maximum length of a string
  117. int max_len_in_all;                                                     //maximum length of a substring which is a substring of every target string
  118. int max_allowed_substr_length;                          //in solution set
  119. int min_allowed_substr_length;                          //in solution set
  120. char toBeCovered[T][M];                                         //target strings
  121. int target_length[T];
  122. double sqrt_m;
  123. ofstream output,output_best;
  124.  
  125. struct node
  126. {
  127. //      int length;
  128.         double cost;
  129.         int substr_idx;
  130.         double pheromone, eta, probability, cum_prob;
  131.         double eta1, eta2;
  132.         void init(int length)
  133.         {
  134.                 pheromone = 10.0;
  135.                 //cost = 1.0;
  136.                 if(length > max_allowed_substr_length)
  137.                         cost = inf;
  138.                 //else if(length > m/4)
  139.                 //      cost = 1.2;
  140.                 else if(length < min_allowed_substr_length)
  141.                         cost = inf;
  142.  
  143. //              else if(length==2)
  144. //                      cost = 1.5;
  145.                 else
  146.                         cost = 1.0;
  147.                 eta = (0.001 * eta1 + 0.999 * eta2 * sqrt(length / sqrt_m) ) / cost;
  148.                 //eta = (0.0995 * eta1 + 0.9 * (eta2) + 0.0005* (length)/(sqrt_m) ) / cost;
  149.                 //eta = (0.001 * eta1 + 0.999 * (eta2) ) / cost; //instances_100_a4_s20-30_l3-10_L50-100
  150.         //      eta = sub_freq[substr_idx]/(id + 0.0);  // needs change
  151.         }
  152. };
  153.                 //targetidx pos
  154. vector<node> adj[T][M];                                 //vector of edges from [T][M]
  155. node all_node[1000005];                                 //just a list of all nodes according to id
  156.  
  157. struct solution
  158. {
  159.         double cost;                                            //cost or fitness of this solution
  160.         set<int> substr_list;                           //set of ids of the substrings
  161.         vector< pair<int,int> > jumps;          //represents the jumps -> <target_string_idx, pos>
  162.         void init()
  163.         {
  164.                 substr_list.clear();
  165.                 jumps.clear();
  166.         }
  167.         void calc_cost()
  168.         {
  169.                 int i;
  170.                 cost = 0.0;
  171.                 vector<int> now_here(substr_list.begin(),substr_list.end());
  172.                 for(i=0;i<now_here.size();i++)
  173.                         cost += all_node[now_here[i]].cost;
  174.         //      printf("calc cost ok \n");
  175.         }
  176. }global_best, local_best, current;
  177.  
  178.  
  179. void addtoTrie(trie *now,int i,int j,int pos)
  180. {
  181. //      if(pos-j+1>13)
  182. //              return;
  183.         if(!toBeCovered[i][pos])
  184.                 return;
  185.  
  186.         pair<int,int> access_target = calc_target_idx(i);
  187.         int temp_idx;  
  188.         int v = toBeCovered[i][pos];
  189.         if(!now->next[v])
  190.         {
  191.                 now->next[v] = new trie();
  192.                 now->next[v]->val = id++;
  193.                
  194.         //      node temp;
  195.         //      temp.substr_idx = now->next[v]->val;
  196.         //      temp.init();
  197.         //      all_node[temp.substr_idx] = temp;
  198.  
  199.         //      now->next[v]->freq = 1;
  200.                 temp_idx = now->next[v]->val;
  201.                 sub_freq[temp_idx] = 0;
  202.                 sub_freq_target[temp_idx] = 0;
  203.                 memset(sub_is_target[temp_idx],0,sizeof(sub_is_target[temp_idx]));
  204.         }
  205.        
  206.         temp_idx = now->next[v]->val;
  207.         sub_freq[temp_idx]++;
  208.         sub_start[temp_idx] = make_pair(i,j);
  209.         sub_end [temp_idx]  = make_pair(i,pos);
  210.         if( !(sub_is_target[temp_idx][access_target.first] & (1<<access_target.second)) )
  211.                 sub_is_target[temp_idx][access_target.first] |= (1<<access_target.second),
  212.                 sub_freq_target[temp_idx]++;
  213. /*
  214.         if( !(now->next[v]->istarget[ access_target.first ] & (1<<access_target.second)))
  215.                 now->next[v]->istarget[ access_target.first ] |= (1<<access_target.second),
  216.                 now->next[v]->freq_target++;
  217. */
  218.         //sub_freq[now->next[v]->val] = now->next[v]->freq;
  219.        
  220.         node temp;
  221.         temp.substr_idx = now->next[v]->val;
  222.         adj[i][j].push_back(temp);
  223.  
  224.         addtoTrie(now->next[v],i,j,pos+1);
  225. }
  226.  
  227. void buildGraph()
  228. {
  229. //      all_node.clear();
  230.         root = new trie();
  231.         root->val = id++;
  232.         int i,j;
  233.         id = 0;
  234.         for(i=0;i<n;i++)
  235.                 for(j=0;toBeCovered[i][j];j++)
  236.                 {
  237.                         adj[i][j].clear();
  238.                         addtoTrie(root,i,j,j);
  239.                 }
  240. //      printf("build graph ok \n");
  241. }
  242.  
  243. void initialize_pheromone()
  244. {
  245.         int i,j,k;
  246.         for(i=0;i<n;i++)
  247.                 for(j=0;toBeCovered[i][j];j++)
  248.                 {
  249.                         for(k=0;k<adj[i][j].size();k++)
  250.                         {
  251.                                 adj[i][j][k].init(k+1);
  252.                                 all_node[adj[i][j][k].substr_idx] = adj[i][j][k];
  253.                         }
  254.                 }
  255.         //      printf("initialze pheromone ok \n");
  256. }
  257.  
  258. /**********************************************/
  259. // may be useful for L-cover
  260. int mark_now;
  261. int target_is_cover;                            // idx of target string for is_cover
  262. //          target_idx from to
  263. char dp_is_cover[T][M][M];                      // is it possible to cover [from, to] within a finite cost?
  264.  
  265. int is_cover(int from, int to)
  266. {
  267.         if(to-from+1 < min_allowed_substr_length)
  268.         {
  269.                 dp_is_cover[target_is_cover][from][to] = mark_now - 1;
  270.                 return dp_is_cover[target_is_cover][from][to];
  271.         }
  272.         if(to-from+1 <= max_allowed_substr_length)
  273.         {
  274.                 dp_is_cover[target_is_cover][from][to] = mark_now;
  275.                 return dp_is_cover[target_is_cover][from][to];
  276.         }
  277.         if(dp_is_cover[target_is_cover][from][to] >= mark_now - 1)
  278.                 return dp_is_cover[target_is_cover][from][to];
  279.         int i, vv;
  280.         dp_is_cover[target_is_cover][from][to] = mark_now - 1;
  281.        
  282.         for(i = from + min_allowed_substr_length - 1; i < to; i++)
  283.         {
  284.                 if( i - from + 1 > max_allowed_substr_length)
  285.                         continue;
  286.                 vv = is_cover(i+1,to);
  287.                 if( vv == mark_now)
  288.                 {
  289.                         dp_is_cover[target_is_cover][from][to] = mark_now;
  290.                         break;
  291.                 }
  292.         }
  293.        
  294.         return dp_is_cover[target_is_cover][from][to];
  295. }
  296.  
  297. void calc_is_cover(void)
  298. {
  299.         mark_now += 2;
  300.         int i,j,k,ffff;
  301.         for(i = 0; i < n; i++)
  302.         {      
  303.                 target_is_cover = i;
  304.                 for(j = 0; toBeCovered[i][j]; j++)
  305.                         for(k = j; toBeCovered[i][k];k++)
  306.                         {
  307.                 //              dp_is_cover[i][j][k] = -1;
  308.                                 ffff = is_cover(j,k);
  309.                         }
  310.         }
  311. }
  312. /******************************************/
  313.  
  314. void visitTrie(trie *now,int i,int j,int pos)
  315. {
  316. //      if(pos-j+1>13)
  317. //              return;
  318.         if(!toBeCovered[i][pos])
  319.                 return;
  320.  
  321.         int v = toBeCovered[i][pos];
  322.         adj[i][j][pos-j].eta1 = sub_freq[now->next[v]->val]/(double)(id+0.0);
  323.         adj[i][j][pos-j].eta2 = sub_freq_target[now->next[v]->val]/(double)(n+0.0);
  324.        
  325.         if(sub_freq_target[now->next[v]->val] == n)
  326.                 max_len_in_all = MAX(max_len_in_all,pos-j+1);
  327.         visitTrie(now->next[v],i,j,pos+1);
  328. }
  329.  
  330. void calc_eta()
  331. {
  332.         max_len_in_all = 0;
  333.         int i, j;
  334.         for(i=0;i<n;i++)
  335.                 for(j=0;toBeCovered[i][j];j++)
  336.                         visitTrie(root,i,j,j);
  337. }
  338.  
  339. int run, ant;
  340.  
  341. int select_jump(int target_idx,int pos,int flag,int interval_count)
  342. {
  343. //      printf("select jump enter \n");
  344.         if(interval_count == (L-1))
  345.                 return (target_length[target_idx] - 1);
  346.         if(target_length[target_idx] - pos <= 2*min_allowed_substr_length)
  347.                 return (target_length[target_idx] - 1);
  348.         int i, ret_idx;
  349.         double sum_heu = 0, max_prob = -1;
  350.         for(i = 0; i < adj[target_idx][pos].size(); i++)
  351.         {
  352.                 if(fabs(adj[target_idx][pos][i].cost - inf) < 1e-5)
  353.                         continue;
  354.                 if(toBeCovered[target_idx][pos+i+1] && dp_is_cover[target_idx][pos+i+1][target_length[target_idx]-1]==(mark_now - 1))
  355.                         continue;
  356.                 double now_pheromone = adj[target_idx][pos][i].pheromone - minus_pheromone;
  357.                 double tabu_factor = 1.0 - Tabu_prob * (sub_tabu[ adj[target_idx][pos][i].substr_idx ] == tabu_active);
  358.         /*      if(flag)
  359.                 {
  360.                         now_pheromone = MAX(t_min,now_pheromone);
  361.                         now_pheromone = MIN(t_max,now_pheromone);
  362.                 } */
  363.                 sum_heu += ( pow(now_pheromone, alpha) * pow(adj[target_idx][pos][i].eta, beta) * tabu_factor);
  364.         }
  365.         for(i = min_allowed_substr_length - 1; i < adj[target_idx][pos].size(); i++)
  366.         {
  367.                 if(fabs(adj[target_idx][pos][i].cost - inf) < 1e-5)
  368.                 {
  369.                         adj[target_idx][pos][i].cum_prob = adj[target_idx][pos][i-1].cum_prob;
  370.                         continue;
  371.                 }
  372.                 if(toBeCovered[target_idx][pos+i+1] && dp_is_cover[target_idx][pos+i+1][target_length[target_idx]-1]==(mark_now - 1))
  373.                 {
  374.                         adj[target_idx][pos][i].cum_prob = adj[target_idx][pos][i-1].cum_prob;
  375.                         continue;
  376.                 }
  377.                 double now_pheromone = adj[target_idx][pos][i].pheromone - minus_pheromone;
  378.                 double tabu_factor = 1.0 - Tabu_prob * (sub_tabu[ adj[target_idx][pos][i].substr_idx ] == tabu_active);
  379.         /*      if(flag)
  380.                 {
  381.                         now_pheromone = MAX(t_min,now_pheromone);
  382.                         now_pheromone = MIN(t_max,now_pheromone);
  383.                 }*/
  384.                 double neumerator = pow(now_pheromone, alpha) * pow(adj[target_idx][pos][i].eta, beta) * tabu_factor;
  385.                 adj[target_idx][pos][i].probability =  neumerator / sum_heu;
  386.  
  387.                 if(adj[target_idx][pos][i].probability > max_prob)
  388.                         max_prob = adj[target_idx][pos][i].probability,
  389.                         ret_idx = i + pos;
  390.  
  391.                 if(i == min_allowed_substr_length - 1)
  392.                         adj[target_idx][pos][i].cum_prob = adj[target_idx][pos][i].probability;
  393.                 else
  394.                         adj[target_idx][pos][i].cum_prob = adj[target_idx][pos][i].probability + adj[target_idx][pos][i-1].cum_prob;
  395.         }
  396.         //int seed;
  397.         //seed = (long int) time(NULL);
  398.  
  399.         if(flag%2 == 0)
  400.         {
  401.                 double q = ran01(&seed);
  402.                 if(q < q0)
  403.                         return ret_idx;
  404.         }
  405.  
  406.         double random_select = ran01(&seed);
  407.  
  408.         for(i=0;i<adj[target_idx][pos].size();i++)
  409.         {
  410.                 if(fabs(adj[target_idx][pos][i].cost - inf) < 1e-5)
  411.                         continue;
  412.                 if(toBeCovered[target_idx][pos+i+1] && dp_is_cover[target_idx][pos+i+1][target_length[target_idx]-1]==(mark_now - 1))
  413.                         continue;
  414.  
  415.                 if(adj[target_idx][pos][i].cum_prob >= (random_select - 1e-12) )
  416.                         return i+pos;
  417.         }
  418. }
  419.  
  420. int select_backjump(int target_idx,int pos,int flag,int interval_count)
  421. {
  422.         if(interval_count == (L-1))
  423.                 return 0;
  424.         if(pos + 1 <= 2*min_allowed_substr_length)
  425.                 return 0;
  426.         int i, ret_idx;
  427.         double sum_heu = 0, max_prob = -1;
  428.         for(i = pos; i >= 0 ; i--)
  429.         {
  430.                 if(fabs(adj[target_idx][i][pos - i].cost - inf) < 1e-5)
  431.                         continue;
  432.                 if(i && dp_is_cover[target_idx][0][i-1]==(mark_now - 1))
  433.                         continue;
  434.                 double now_pheromone = adj[target_idx][i][pos - i].pheromone - minus_pheromone;
  435.                 double tabu_factor = 1.0 - Tabu_prob * (sub_tabu[ adj[target_idx][i][pos - i].substr_idx ] == tabu_active);
  436.        
  437.                 sum_heu += ( pow(now_pheromone, alpha) * pow(adj[target_idx][i][pos - i].eta, beta) * tabu_factor);
  438.         }
  439.         for(i = pos; i >= 0 ; i--)
  440.         {
  441.                 if(fabs(adj[target_idx][i][pos - i].cost - inf) < 1e-5)
  442.                         continue;
  443.                 if(i && dp_is_cover[target_idx][0][i-1]==(mark_now - 1))
  444.                         continue;
  445.                 double now_pheromone = adj[target_idx][i][pos - i].pheromone - minus_pheromone;
  446.                 double tabu_factor = 1.0 - Tabu_prob * (sub_tabu[ adj[target_idx][i][pos - i].substr_idx ] == tabu_active);
  447.                 double numerator = ( pow(now_pheromone, alpha) * pow(adj[target_idx][i][pos - i].eta, beta) * tabu_factor);
  448.                 adj[target_idx][i][pos - i].probability = numerator/sum_heu;
  449.  
  450.                 if(adj[target_idx][i][pos - i].probability > max_prob)
  451.                 {
  452.                         max_prob = adj[target_idx][i][pos - i].probability;
  453.                         ret_idx = i;
  454.                 }
  455.  
  456.                 if(pos - i + 1 == min_allowed_substr_length)
  457.                         adj[target_idx][i][pos - i].cum_prob = adj[target_idx][i][pos - i].probability;
  458.                 else
  459.                         adj[target_idx][i][pos - i].cum_prob = adj[target_idx][i+1][pos - i - 1].cum_prob + adj[target_idx][i][pos - i].probability;
  460.         }
  461.        
  462.         if(flag%2 == 0)
  463.         {
  464.                 double q = ran01(&seed);
  465.                 if(q < q0)
  466.                         return ret_idx;
  467.         }
  468.  
  469.         double random_select = ran01(&seed);
  470.         for(i = pos; i >= 0 ; i--)
  471.         {
  472.                 if(fabs(adj[target_idx][i][pos - i].cost - inf) < 1e-5)
  473.                         continue;
  474.                 if(i && dp_is_cover[target_idx][0][i-1]==(mark_now - 1))
  475.                         continue;
  476.                 if(adj[target_idx][i][pos - i].cum_prob >= random_select - 1e-12)
  477.                         return i;
  478.         }
  479.  
  480.         return ret_idx;
  481. }
  482.  
  483. int construction_failure;                       //not used
  484.  
  485. void construction(int ant, int flag)
  486. {
  487. //      printf("construction enter \n");
  488.         //int seed;
  489.         //seed = (int) time(NULL);
  490.         int start_target = (int)(ran01(&seed)*n);
  491.         start_target %= n;
  492.         int i, pos;
  493.         current.init();
  494.         for(i = start_target;;)
  495.         {
  496.                 if(rand()%2)
  497.                 {
  498.                         current.jumps.push_back(make_pair(i,0));
  499.                         pos = 0;
  500.                         int interval_count = 0;                                 //can not be greater than L for L-cover
  501.                         while(toBeCovered[i][pos])
  502.                         {
  503.                                 int k = select_jump(i,pos,flag,interval_count);
  504.                 //              printf("select jump ok \n");
  505.                                
  506.                                 current.substr_list.insert(adj[i][pos][k-pos].substr_idx);
  507.                                 current.jumps.push_back(make_pair(i,k+1));
  508.                                 pos = k+1;
  509.                                 interval_count++;
  510.                         }
  511.                 }
  512.                 else
  513.                 {
  514.                         vector< pair<int,int> > tempv;
  515.                         pos = target_length[i] - 1;
  516.                         tempv.push_back(make_pair(i,pos+1));
  517.                         int interval_count = 0;                                 //can not be greater than L for L-cover
  518.                         while(pos>=0)
  519.                         {
  520.                                 int k = select_backjump(i,pos,flag,interval_count);
  521.                                 current.substr_list.insert(adj[i][k][pos-k].substr_idx);
  522.                                 tempv.push_back(make_pair(i,k));
  523.                                 pos = k-1;
  524.                                 interval_count++;
  525.                         }
  526.  
  527.                         reverse(tempv.begin(),tempv.end());
  528.                         for(int k = 0; k < tempv.size(); k++)
  529.                                 current.jumps.push_back(tempv[k]);
  530.                 }
  531.                 i = (i+1)%n;
  532.                 if(i==start_target)
  533.                         break;
  534.         }
  535.         current.calc_cost();
  536. //      printf("construntion ok \n");
  537. }
  538.  
  539. void backward_construction(int ant, int flag)
  540. {
  541.        
  542. }
  543.  
  544. void update_pheromone(int flag, solution got)                   //need to modify according to flag
  545. {
  546.         //      printf("update pheromone enter\n");
  547.         int i, j, k;
  548.         double k1;
  549.         //      fixing t_max and t_min for the solution we got here
  550.         double tmax_d = evaporation_rate * got.cost;
  551.         t_max = 1.0/tmax_d;
  552.         k1 = m/2.0;                                     // max length of string/2 -> average length jump
  553.         double p_dec = pow(p_best,(double) 1.0/(double)m);
  554.         double tmin_d = (k1-1)*p_dec;
  555.         t_min = t_max*(1-p_dec)/tmin_d;
  556.  
  557.         //      evaporate
  558.         for(i=0;i<n;i++)
  559.                 for(j=0;toBeCovered[i][j];j++)
  560.                         for(k = 0;k<adj[i][j].size();k++)
  561.                         {
  562.                                 adj[i][j][k].pheromone *= (1 - evaporation_rate);
  563.                                 adj[i][j][k].pheromone = MAX(adj[i][j][k].pheromone, t_min);
  564.                                 adj[i][j][k].pheromone = MIN(adj[i][j][k].pheromone, t_max);
  565.                         }
  566.        
  567.         // increase pheromone
  568.         for(i=0;i<got.jumps.size();i++)
  569.         {
  570.                 int uuu = got.jumps[i].first;                   //target string index
  571.                 int vvv = got.jumps[i].second;                  //current position
  572.                 if(!toBeCovered[ uuu ][ vvv ])                  //reached end of this target string
  573.                         continue;      
  574.                 if(i+1>=got.jumps.size())
  575.                         continue;
  576.                 j = got.jumps[i+1].second - 1;
  577.                  //jump position from current position
  578.                 //if(tabu_done)
  579.                 //{
  580.                 //      adj[ uuu ][ vvv ][ j - vvv].pheromone -= evaporation_rate * (1.0/got.cost);
  581.                 //}
  582.                 //else
  583.                         adj[ uuu ][ vvv ][ j - vvv].pheromone += evaporation_rate * (1.0/got.cost);
  584.  
  585.                 adj[ uuu ][ vvv ][ j - vvv].pheromone = MIN(adj[ uuu ][ vvv ][ j - vvv].pheromone, t_max);
  586.                 adj[ uuu ][ vvv ][ j - vvv].pheromone = MAX(adj[ uuu ][ vvv ][ j - vvv].pheromone, t_min);
  587.         }
  588.         //      printf("update pheromone ok \n");
  589. }
  590.  
  591. void output_strings_best(int run)
  592. {
  593.         output_best<<"Run "<<run<<endl;
  594.         int i;
  595. /*      for(i=0;i<global_best.jumps.size();i++)
  596.         {
  597.                 int now_taget_idx = global_best.jumps[i].first;
  598.                 if(i+1>=global_best.jumps.size() || global_best.jumps[i+1].first!=now_taget_idx)
  599.                 {
  600.                         output_best<<endl;
  601.                         continue;
  602.                 }
  603.                 for(int j=global_best.jumps[i].second;j<global_best.jumps[i+1].second;j++)
  604.                         output_best<<toBeCovered[now_taget_idx][j];
  605.                 output_best<<" ";
  606.         }*/
  607.         output_best<<"Solution Set"<<endl;
  608.         vector<int> idx_list(global_best.substr_list.begin(),global_best.substr_list.end());
  609.         for(i=0;i<idx_list.size();i++)
  610.         {
  611.                 int now_idx = idx_list[i];
  612.                 int target_now = sub_start[now_idx].first;
  613.                 for(int j = sub_start[now_idx].second; j<= sub_end[now_idx].second; j++)
  614.                         output_best<<toBeCovered[target_now][j];
  615.                 output_best<<endl;
  616.         }
  617. }
  618.  
  619. solution local_search(solution got)
  620. {
  621.         int tweak = n;
  622.         solution best_local;
  623.         best_local.cost = 1000000000;
  624.         while(tweak--)
  625.         {
  626.                 solution temp = got;
  627.                 int no_remove = 10;
  628.                 while(no_remove--)
  629.                 {
  630.                         int target = (int)(ran01(&seed)*n);
  631.                         target %= n;
  632.                         int i = 0;
  633.                         while(temp.jumps[i].first != target)
  634.                                 i++;
  635.                         int no_partition = 0, j = i+1;
  636.                         while(j < temp.jumps.size() && temp.jumps[j].first==target)
  637.                                 no_partition++,j++;
  638.  
  639.                         int remove_idx = (int)(ran01(&seed)*no_partition);
  640.                         remove_idx %= no_partition;
  641.  
  642.                         if(remove_idx == 0 || remove_idx == no_partition)
  643.                         {
  644.                                 no_remove++;
  645.                                 continue;
  646.                         }
  647.  
  648.                         temp.jumps.erase(temp.jumps.begin()+i+remove_idx);
  649.  
  650.                         temp.substr_list.clear();
  651.  
  652.                         for(i=0;i<temp.jumps.size();i++)
  653.                         {
  654.                                 if(i+1>=temp.jumps.size())
  655.                                         continue;
  656.                                 if(!toBeCovered[temp.jumps[i].first][temp.jumps[i].second])
  657.                                         continue;
  658.                                 j = temp.jumps[i+1].second - 1 - temp.jumps[i].second;
  659.                                 temp.substr_list.insert(adj[temp.jumps[i].first][temp.jumps[i].second][j].substr_idx);
  660.                         }
  661.  
  662.                         temp.calc_cost();
  663.  
  664.                         if(temp.cost<best_local.cost)
  665.                                 best_local = temp;
  666.                 }
  667.         }
  668.         return best_local;
  669. }
  670.  
  671. /*
  672. // for sorting strings according to size, then alphabetical order
  673. bool comp(string founda,string foundb)
  674. {
  675.         if(founda.size() == foundb.size())
  676.                 return founda < foundb;
  677.         return founda.size() < foundb.size();
  678. }
  679. */
  680.  
  681. // for sorting idx_list according to eta
  682. bool comp_substring_eta(int id1, int id2)
  683. {
  684.         return all_node[id1].eta < all_node[id2].eta;
  685. }
  686.  
  687. void mark_tabu2(void)
  688. {
  689.         vector<int> idx_list(global_best.substr_list.begin(),global_best.substr_list.end());
  690.         sort(idx_list.begin(), idx_list.end(), comp_substring_eta);
  691.         int to_be_tabu = idx_list.size()/3, i;
  692.         for(i = 0; i < to_be_tabu; i++)
  693.                 sub_tabu[ idx_list[i] ] = tabu_active;
  694. }
  695.  
  696. bool comp_substring_pheromone(int id1, int id2)
  697. {
  698.         return all_node[id1].pheromone > all_node[id2].pheromone;
  699. }
  700.  
  701. void mark_tabu3(void)
  702. {
  703.         vector<int> idx_list(global_best.substr_list.begin(),global_best.substr_list.end());
  704.         sort(idx_list.begin(), idx_list.end(), comp_substring_pheromone);
  705.         int to_be_tabu = idx_list.size(), i;
  706.         for(i = 0; i < to_be_tabu; i++)
  707.                 sub_tabu[ idx_list[i] ] = tabu_active;
  708. }
  709.  
  710. void mark_tabu(void)
  711. {
  712.         vector<int> idx_list(global_best.substr_list.begin(),global_best.substr_list.end());
  713.         vector<string> string_list;
  714.         int i,j;
  715.         for(i = 0; i < idx_list.size(); i++)
  716.         {
  717.                 int temp_idx = idx_list[i];
  718.                 int target_idx = sub_start[temp_idx].first;
  719.                 char temp_s[1000];
  720.                 for(j = sub_start[temp_idx].second; j <= sub_end[temp_idx].second; j++)
  721.                         temp_s[ j - sub_start[temp_idx].second ] = toBeCovered[target_idx][j];
  722.                 temp_s[ j - sub_start[temp_idx].second ] = 0;
  723.  
  724.                 string_list.push_back(temp_s);
  725.         }
  726. //      sort korle accordingly idx_list o update korte hobe.   
  727. //      sort(string_list.begin(),string_list.end(),comp());
  728.        
  729.         int cnt_tabu = 0;
  730.         for(i = 0; i < string_list.size(); i++)
  731.                 for(j = 0; j < string_list.size(); j++)
  732.                 {
  733.                         if(j == i)
  734.                                 continue;
  735.                         if(string_list[i].size() == 1)
  736.                                 continue;
  737.                         int vvv = (int)string_list[j].find(string_list[i]);
  738.                         if(vvv != -1)                                                                                           // found i in j
  739.                         {
  740.                 //              if(rand()%2)
  741.                 //                      continue;
  742.                                 if(rand()%5)
  743.                                         sub_tabu[ idx_list[i] ] = tabu_active;                                  // tabu it
  744.                                 else
  745.                                         sub_tabu[ idx_list[j] ] = tabu_active;
  746.                                 cnt_tabu++;
  747.                         }                      
  748.                 }
  749.  
  750.                 if(cnt_tabu <= idx_list.size()/5)
  751.                         mark_tabu2();
  752. }
  753.  
  754. void ACO(int flag)                                      // flag - 0 - ACS, 1 - MMAS, 2 - Hybrid
  755. {
  756.         calc_eta();
  757.         initialize_pheromone();
  758.        
  759.         calc_is_cover();       
  760.         int iter,iter_best;
  761.         double best_time;
  762.        
  763.         for(run = 1; run <= MAX_RUN; run++)
  764.         {
  765.         /*      minus_pheromone = 0.0;  -> not used */
  766.                 tabu_done = 0;
  767.                 global_best.cost = 100000000;
  768.                 double st_time = time(NULL);
  769.                 int idx_fgb = 0;
  770.                 initialize_pheromone();
  771.                 iter = 0;
  772.                 int no_change = 0;
  773.  
  774.                 int last_tabu_on = 0;
  775.                 int local_best_all = 100000000;
  776.                 while(1)
  777.                 {
  778.                         local_best.cost = 100000000;
  779.                         for(ant = 1; ant <= POPULATION_SIZE; ant++)
  780.                         {
  781.                                 srand(iter*iter + rand()%100);
  782.                                 construction(ant, flag);
  783.                         //      solution from_local = local_search(current);
  784.                         //      if(from_local.cost < current.cost)
  785.                         //              current = from_local;
  786.                 //              printf("%lf\n",current.cost);
  787.                                 if(current.cost<local_best.cost)
  788.                                         local_best = current;
  789.                 //              if(flag == 0)                                           // for ACS
  790.                 //                      update_pheromone(flag);
  791.                         }
  792.                         if(local_best.cost == local_best_all)
  793.                                 no_change++;
  794.                         else{
  795.                                 local_best_all = local_best.cost;
  796.                                 no_change = 0;
  797.                         }
  798. /*
  799.                         if(rand()%7 == 0)
  800.                         {
  801.                                 solution from_local = local_search(local_best);
  802.                                 if(from_local.cost < local_best.cost)
  803.                                         local_best = from_local;
  804.                         }
  805.                         */
  806.                         if(local_best.cost < global_best.cost)
  807.                         {
  808.                                 //solution from_local = local_search(local_best);
  809.                                 //if(from_local.cost < local_best.cost)
  810.                                 //      local_best = from_local;
  811.                                 global_best = local_best;
  812.                                 double nd_time = time(NULL);
  813.                                 double elapsed = (nd_time - st_time);
  814.                                 iter_best = iter;
  815.                                 best_time = elapsed;
  816.                                 //no_change = 0;
  817.                         }
  818.                         //else
  819.                         //      no_change++;
  820.  
  821.                         //schedule for tabu
  822.                         if(tabu_done > 0)
  823.                                 tabu_done++;
  824.                        
  825.                 /*      if((no_change > NO_UPDATE/4) && (tabu_done==0) && ((iter - last_tabu_on) > NO_UPDATE/4))
  826.                         {
  827.                                 printf("here tabu on iteration %d\n tabu_idx %d\n",iter,tabu_active);
  828.                                 tabu_done = 1;
  829.                                
  830.                                 if(rand()%2 == 0)
  831.                                         mark_tabu();
  832.                                 //else
  833.                                 {
  834.                                         if(rand()%3 == 0)
  835.                                                 mark_tabu2();
  836.                                         //else
  837.                                                 mark_tabu3();
  838.                                 }
  839.                         }
  840.  
  841.                         if(tabu_done > NO_UPDATE/2)
  842.                         {
  843.                                 printf("here tabu off %d iteration %d\n",tabu_done,iter);
  844.                                 tabu_done = 0;
  845.                        
  846.                                 tabu_active++;
  847.                         }
  848.                        
  849.                         if(tabu_done)
  850.                                 last_tabu_on = iter;*/
  851.                         // FOR both ACS and MMAS
  852.                         // Scheduling
  853.                         int fgb[] = {7,6,5,4,3,2,1};
  854.                         int iter_fgb[] = {25,50,100,200,300,400,10000000};
  855.                        
  856.                 /*      if(iter > 120 && iter%50 == 0)
  857.                         {
  858.                                 beta -= alpha;
  859.                                 if(beta < alpha)
  860.                                         beta = alpha + 1;
  861.                         }*/
  862.                         if(iter > iter_fgb[idx_fgb])
  863.                                 idx_fgb++;
  864.  
  865.                         if(iter%fgb[idx_fgb]==0 && iter>150)
  866.                                 update_pheromone(flag, global_best);                        
  867.                         else
  868.                                 update_pheromone(flag, local_best);
  869.                
  870.                         iter++;
  871.                        
  872.                         printf("iteration %d global_best %lf %d local_best %lf\n",iter,global_best.cost,global_best.substr_list.size(),local_best.cost);
  873.                         if(no_change > NO_UPDATE)
  874.                         {
  875.                                 // record
  876.                                 break;
  877.                         }
  878.  
  879.                         double nd_time = time(NULL);
  880.                         double elapsed = (nd_time - st_time);
  881.                         if(elapsed > allowed_time)
  882.                         {
  883.                                 // record
  884.                                 break;
  885.                         }
  886.                 }
  887.                 output<<"Run "<<run<<endl;
  888.                 output<<iter_best<<" "<<best_time<<" "<<global_best.substr_list.size()<<endl;
  889.                 output_strings_best(run);
  890.         }
  891. }
  892.  
  893. int readInput(int fileidx)
  894. {
  895.         /*if(scanf("%d",&n)!=1)                                         // number of target strings
  896.                 return 0;*/
  897.         int i;
  898.         m = 0;
  899.         n = 10; //all set length 100
  900.         for(i=0;i<n;i++)
  901.         {
  902.                 scanf("%s",toBeCovered[i]);
  903.                 int v = strlen(toBeCovered[i]);
  904.                 m = MAX(v,m);
  905.                 target_length[i] = v;
  906.         }
  907.         sqrt_m = sqrt((double)m);
  908.         max_allowed_substr_length = m/2;
  909.         min_allowed_substr_length = 3;
  910.         return 1;
  911. }
  912.  
  913. void parameterRead(int fileidx)
  914. {
  915.         char str[55];
  916.         sprintf(str,"input\\param%.2d.txt",fileidx);
  917.         freopen(str,"r",stdin);
  918.         scanf("%lf%lf%lf",&alpha,&beta,&evaporation_rate);
  919.         scanf("%d",&L);
  920.         scanf("%lf",&allowed_time);
  921. //      fclose(fp);
  922. }
  923.  
  924. void init_output_file(int fileidx)
  925. {
  926.         char str[200];
  927.         sprintf(str,"output\\%s\\outputtest%.2d.txt",str_folder,fileidx);
  928.         output.open(str);
  929.         sprintf(str,"output\\%s\\outputbest%.2d.txt",str_folder,fileidx);
  930.         output_best.open(str);
  931. }
  932.  
  933. int main()
  934. {
  935.         int fileidx;
  936.         //scanf("%d",&fileidx);
  937.         char str_command[100];
  938.         strcpy(str_folder,"instances_10_a50_s2-10_l3-10_L50-100");
  939.         sprintf(str_command,"mkdir output\\%s",str_folder);
  940.         system(str_command);
  941.         for(fileidx = 2; fileidx <= 50; fileidx+=2)
  942.         {
  943.                 //parameterRead(fileidx);
  944.                        
  945.                 init_output_file(fileidx);
  946.                 char str[200];
  947.                
  948.                 sprintf(str,"instances\\%s\\instance%.2d",str_folder,fileidx);
  949.                 freopen(str,"r",stdin);
  950.         //      while(readInput(fileidx))
  951.                 {
  952.                         readInput(fileidx);
  953.                         buildGraph();
  954.                         alpha = 5.0;
  955.                         beta = 20.0;
  956.                         evaporation_rate = 0.02;
  957.                         L = 200;
  958.                         allowed_time = 7200;
  959.  
  960.                         output<<"alpha "<<alpha<<" beta "<<beta<<" evaporation_rate "<<evaporation_rate<<" with tabu"<<endl;
  961.                         output_best<<"alpha "<<alpha<<" beta "<<beta<<" evaporation_rate "<<evaporation_rate<<" with tabu"<<endl;
  962.                         tabu_active++;
  963.                         ACO(2);
  964.         /*              beta += 0.5;
  965.                         output<<"alpha "<<alpha<<"beta "<<beta<<"evaporation_rate "<<evaporation_rate<<endl;
  966.                         output_best<<"alpha "<<alpha<<"beta "<<beta<<"evaporation_rate "<<evaporation_rate<<endl;
  967.                         ACO(2);
  968.                         evaporation_rate -= 0.005;
  969.                         output<<"alpha "<<alpha<<"beta "<<beta<<"evaporation_rate "<<evaporation_rate<<endl;
  970.                         output_best<<"alpha "<<alpha<<"beta "<<beta<<"evaporation_rate "<<evaporation_rate<<endl;
  971.                         ACO(2);*/
  972.                 }
  973.                 output.close();
  974.                 output_best.close();
  975.         }
  976.         return 0;
  977. }
clone this paste RAW Paste Data