fork download
copy 
  1. #include <stdio.h>
  2. #include <malloc.h>
  3. #include <math.h>
  4. #include <limits.h>
  5. #include <stdlib.h>
  6. #include <time.h>
  7.  
  8. struct Set
  9. {
  10.   int n;    /* number of elements in the set */
  11.   int min;  /* minimum distance between elements of this set */
  12.   int best; /* best minimum distance for all sets up to this set */
  13. };
  14.  
  15. /* Set sizes */
  16. struct Set sets[] = {
  17.   { 19 },
  18.   { 19 },
  19.   { 18 },
  20.   { 16 },
  21.   { 16 },
  22.   { 14 },
  23.   { 11 },
  24.   { 11 },
  25.   { 11 },
  26.   { 9 },
  27.   { 9 },
  28.   { 7 },
  29.   { 7 },
  30.   { 7 },
  31.   { 6 },
  32.   { 5 },
  33.   { 4 },
  34.   { 3 },
  35.   { 2 },
  36.   { 2 },
  37. };
  38. #define Singles 6 /* number of one-element sets */
  39.  
  40. #define RT_BEGIN 10.        /* starting reciprocal "temperature" */
  41. #define RT_END 1000.        /* ending reciprocal "temperature" */
  42. #define RT_STEP 3.0e-5      /* speed of "temperature" change */
  43. #define RT_FINE_STEP 1.0e-5 /* speed of "temperature" change on final step*/
  44. #define IMPROVE_LIMIT 0.5/*portion of array values with harder qality function*/
  45. #define SRAND 1             /* randomize */
  46. #define SKIP_GROUPS 0       /* perform simulated annealing only once */
  47.  
  48. #define M (sizeof(sets) / sizeof(struct Set))
  49. #define SinglesSet (-1)
  50. #define INF INT_MAX
  51.  
  52. #define SWAP(a, b) { struct A t = a; a = b; b = t; }
  53.  
  54. struct A
  55. {
  56.   int set; /* set number */
  57.   int l;   /* distance to the nearest element of the same set to the left */
  58.   int r;   /* distance to the nearest element of the same set to the right */
  59. };
  60.  
  61. struct Align
  62. {
  63.   struct A* a;    /* array */
  64.   double* l;      /* pre-computed table of logarithms */
  65.   int* m;         /* histograms of distances */
  66.   int n;          /* number of elements in the array */
  67.   int nn;         /* number of elements in largest sets */
  68.   int bestMin;    /* optimal min.distance for largest sets */
  69.   int startVar;   /* number of largest sets */
  70.   int startWatch; /* start of not-yet optimized portion of sets */
  71. };
  72.  
  73. /* Quality function */
  74. static double eval(double* l, int n, int x)
  75. {
  76.   ++x;
  77.  
  78.   if (x > 0)
  79.   {
  80.     return (x >= n)? 0: l[x];
  81.   }
  82.   else
  83.   {
  84.     return -(1 << (1 - x));
  85.   }
  86. }
  87.  
  88. /* Calculates how much quality may be improved after swapping elements
  89.    left and left+1 */
  90. static double compare(struct Align* a, int left)
  91. {
  92.   double result = 0.;
  93.   int right = left + 1;
  94.  
  95.   if (a->a[left].set != SinglesSet)
  96.   {
  97.     int best = sets[a->a[left].set].best;
  98.  
  99.     result +=
  100.       eval(a->l, a->n, a->a[left].l + 1 - best) -
  101.       eval(a->l, a->n, a->a[left].l - best) +
  102.       eval(a->l, a->n, a->a[left].r - 1 - best) -
  103.       eval(a->l, a->n, a->a[left].r - best);
  104.   }
  105.  
  106.   if (a->a[right].set != SinglesSet)
  107.   {
  108.     int best = sets[a->a[right].set].best;
  109.  
  110.     result +=
  111.       eval(a->l, a->n, a->a[right].l - 1 - best) -
  112.       eval(a->l, a->n, a->a[right].l - best) +
  113.       eval(a->l, a->n, a->a[right].r + 1 - best) -
  114.       eval(a->l, a->n, a->a[right].r - best);
  115.   }
  116.  
  117.   return result;
  118. }
  119.  
  120. /* Harden quality function for not-yet-optimized portion of sets if possible */
  121. static void improve(struct Align* a)
  122. {
  123.   int i;
  124.   for (i = 0; i < a->startWatch; ++i)
  125.   {
  126.     if (sets[i].min < sets[i].best)
  127.       return;
  128.   }
  129.  
  130.   for (i = a->startWatch; i < M; ++i)
  131.   {
  132.     if (sets[i].min <= sets[i].best)
  133.       return;
  134.   }
  135.  
  136.   for (i = a->startWatch; i < M; ++i)
  137.   {
  138.     ++sets[i].best;
  139.   }
  140. }
  141.  
  142. /* Update distances for one of the swapped elements,
  143.   update minimal distance for correspondent set */
  144. static void
  145. update(struct Align* a, int* d, int theSet, int pos, int diff, int lr)
  146. {
  147.   int better = 0;
  148.   --a->m[a->n * theSet + *d];
  149.  
  150.   if (sets[theSet].min == *d &&
  151.       a->m[a->n * theSet + *d] == 0)
  152.   {
  153.     sets[theSet].min += diff;
  154.     better = (diff == 1)? 1: 0;
  155.   }
  156.  
  157.   *d += diff;
  158.   pos += *d * diff * lr * -1;
  159.   ++a->m[a->n * theSet + *d];
  160.  
  161.   if (*d == sets[theSet].min - 1)
  162.   {
  163.     sets[theSet].min = *d;
  164.   }
  165.  
  166.   if (diff * lr == 1) /* left */
  167.   {
  168.     a->a[pos].r += diff;
  169.   }
  170.   else /* right */
  171.   {
  172.     a->a[pos].l += diff;
  173.   }
  174.  
  175.   if (better)
  176.   {
  177.     improve(a);
  178.   }
  179. }
  180.  
  181. /* Swap elements left and left+1 */
  182. static void swap(struct Align* a, int left)
  183. {
  184.   int right = left + 1;
  185.   SWAP(a->a[left], a->a[right])
  186.  
  187.   if (a->a[left].l < a->n)
  188.     update(a, &(a->a[left].l), a->a[left].set, left, -1, -1);
  189.  
  190.   if (a->a[right].r < a->n)
  191.     update(a, &(a->a[right].r), a->a[right].set, right, -1, 1);
  192.  
  193.   if (a->a[left].r < a->n)
  194.     update(a, &(a->a[left].r), a->a[left].set, left, 1, -1);
  195.  
  196.   if (a->a[right].l < a->n)
  197.     update(a, &(a->a[right].l), a->a[right].set, right, 1, 1);
  198. }
  199.  
  200. /* Main loop of simulated annealing */
  201. static void iterate(struct Align* a, double step)
  202. {
  203.   int i;
  204.   double rt = RT_BEGIN;
  205.  
  206.   while (rt < RT_END)
  207.   {
  208.     /*printf("\nrt=%g\n", rt);*/
  209.  
  210.     for (i = 0; i != a->n; ++i)
  211.     {
  212.       double d;
  213.       int left;
  214.  
  215.       do {
  216. 	left  = rand() % (a->n - 1);
  217.       } while (a->a[left].set == a->a[left+1].set);
  218.  
  219.       d = compare(a, left);
  220.  
  221.       if (d >= 0. || exp(d * rt) * RAND_MAX > rand())
  222.       {
  223. 	swap(a, left);
  224.       }
  225.     }
  226.  
  227.     rt *= 1.0 + step;
  228.   }
  229. }
  230.  
  231. /* Perform simulated annealing for equally-sized groups of sets,
  232.   which allows to incrementally tighten quality function for these groups */
  233. static void groups(struct Align* a)
  234. {
  235.   for (a->startWatch = a->startVar; a->startWatch < M;)
  236.   {
  237.     int oldN = sets[a->startWatch].n;
  238.     iterate(a, RT_STEP);
  239.  
  240. #if SKIP_GROUPS
  241.     return;
  242. #endif
  243.  
  244.     for (; a->startWatch < M && sets[a->startWatch].n == oldN; ++a->startWatch)
  245.     {
  246.       a->nn += sets[a->startWatch].n;
  247.  
  248.       if (a->nn > a->n * IMPROVE_LIMIT)
  249.       {
  250. 	goto exit;
  251.       }
  252.     }
  253.   }
  254.  
  255.  exit:
  256.   iterate(a, RT_FINE_STEP);
  257. }
  258.  
  259. /* Allocate memory, initialize everything */
  260. static void setup()
  261. {
  262.   int i;
  263.   int theSet = 0;
  264.   int theSetElems = sets[theSet].n;
  265.   struct Align a;
  266.  
  267.   a.a = 0;
  268.   a.l = 0;
  269.   a.m = 0;
  270.   a.n = Singles;
  271.   a.startVar = 0;
  272.  
  273.   for (i = 0; i < M; ++i)
  274.   {
  275.     a.n += sets[i].n;
  276.  
  277.     if (!a.startVar && sets[i].n != sets[0].n)
  278.     {
  279.       a.startVar = i;
  280.     }
  281.  
  282.     if (!a.startVar)
  283.     {
  284.       a.nn += sets[i].n;
  285.     }
  286.   }
  287.  
  288.   a.a = calloc(a.n, sizeof(struct A));
  289.   a.l = calloc(a.n, sizeof(double));
  290.   a.m = calloc(a.n, M * sizeof(int));
  291.  
  292.   for (i = 1; i < a.n; ++i)
  293.   {
  294.     a.l[i] = log(i);
  295.   }
  296.  
  297.   a.bestMin = (a.n - a.startVar) / (sets[0].n - 1);
  298.   i = 0;
  299.   while (1)
  300.   {
  301.     a.a[i].set = theSet;
  302.     a.a[i].l = (theSetElems == sets[theSet].n)? INF: 1;
  303.     --theSetElems;
  304.     a.a[i].r = theSetElems? 1: INF;
  305.     ++i;
  306.  
  307.     if (!theSetElems)
  308.     {
  309.       sets[theSet].min = 1;
  310.       sets[theSet].best = a.bestMin;
  311.       a.m[a.n * theSet + 1] = sets[theSet].n - 1;
  312.       ++theSet;
  313.  
  314.       if (theSet == M)
  315. 	break;
  316.  
  317.       theSetElems = sets[theSet].n;
  318.     }
  319.   }
  320.   for (; i < a.n; ++i)
  321.   {
  322.     a.a[i].set = SinglesSet;
  323.     a.a[i].l = INF;
  324.     a.a[i].r = INF;
  325.   }
  326.  
  327.   groups(&a);
  328.  
  329.   /*for (i = 0; i < a.n; ++i)
  330.   {
  331.     printf("%d ", a.a[i].set);
  332.   }*/
  333.  
  334.   free(a.m);
  335.   free(a.l);
  336.   free(a.a);
  337. }
  338.  
  339. int main()
  340. {
  341. #if SRAND
  342.   srand(time(0));
  343. #endif
  344.  
  345.   setup();
  346.  
  347.   /* Print set size, minimal distance, limit for minimal distance */
  348.   int i;
  349.   for (i = 0; i < M; ++i)
  350.   {
  351.     printf("\n%d\t%d", sets[i].n, sets[i].min);
  352.  
  353.     if (sets[i].best != sets[i].min)
  354.     {
  355.       printf("\t%d", sets[i].best);
  356.     }
  357.   }
  358.   puts("\n");
  359.  
  360.   return 0;
  361. }
  362.  
Not running #stdin #stdout 0s 0KB
comments (?)
stdin
copy 
Standard input is empty
stdout
copy 
Standard output is empty
https://ideone.com/OGkkc
language:
C (gcc 8.3)
created:
11 years ago
visibility:
secret

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!