min_dnf.c

// gcc min_dnf.c -o min_dnf -lgmp

#include <gmp.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#define MAX (1 << 21)

typedef struct PAIR {
  int fst;
  int snd;
} pair;

typedef struct PNODE {
  pair val;
  struct PNODE* prev;
  struct PNODE* next;
} pnode;

typedef struct PLIST {
  int size;
  pnode* head;
  pnode* tail;
} plist;

void push_front(plist* l, pair x) {
  pnode* tmp = (pnode*)calloc(1, sizeof(pnode));
  tmp->val.fst = x.fst;
  tmp->val.snd = x.snd;
  tmp->prev = NULL;
  tmp->next = NULL;
  if (l->head != NULL) {
    tmp->next = l->head;
    l->head->prev = tmp;
  } else {
    l->tail = tmp;
  }
  l->head = tmp;
  l->size += 1;
}

/* void push_back(plist* l, pair x) { */
/*   pnode* tmp = (pnode*)calloc(1, sizeof(pnode)); */
/*   tmp->val.fst = x.fst; */
/*   tmp->val.snd = x.snd; */
/*   if (l->head != NULL) { */
/*     tmp->next = l->head; */
/*   } */
/*   l->head = tmp; */
/*   l->size += 1; */
/* } */

void push_front_u(plist* l, pair x) {
  for (pnode* p = l->head; p != NULL; p = p->next) {
    if (x.fst == p->val.fst && x.snd == p->val.snd) {
      return;
    }
  }
  
  pnode* tmp = (pnode*)calloc(1, sizeof(pnode));
  tmp->val.fst = x.fst;
  tmp->val.snd = x.snd;
  tmp->prev = NULL;
  if (l->head != NULL) {
    tmp->next = l->head;
    l->head->prev = tmp;
  } else {
    l->tail = tmp;
  }
  l->head = tmp;
  l->size += 1;
}

int delete(plist* l, pair x) {
  if (l == NULL || l->size == 0) {
    return -1;
  }
  pnode *p, *q;
  if (l->head->val.fst == x.fst && l->head->val.snd == x.snd) {
    p = l->head;
    l->head = l->head->next;
    free(p);
    l->size -= 1;
    return 0;
  }
  for (p = l->head; p != NULL; ) {
    q = p;
    p = p->next;
    if (p != NULL && p->val.fst == x.fst && p->val.snd == x.snd) {
      q->next = p->next;
      l->size -= 1;
      free(p);
      return 0;
    }
  }
  return -1;
}

void clear(plist* l) {
  if (l == NULL) {
    return;
  }
  pnode *p, *q;
  for (p = l->head; p != NULL; ) {
    q = p;
    p = p->next;
    free(q);
  }
  l->size = 0;
  l->head = NULL;
}

void plist_print(plist* l) {
  if (l == NULL) {
    return;
  }
  for (pnode* p = l->head; p != NULL; p = p->next) {
    printf("(%d, %d) ", p->val.fst, p->val.snd);
  }
}

void plist_print_reverse(plist* l) {
  if (l == NULL) {
    return;
  }
  for (pnode* p = l->tail; p != NULL; p = p->prev) {
    printf("(%d, %d) ", p->val.fst, p->val.snd);
  }
}


int get_number_of_variables(char* tt) {
  long unsigned x = strlen(tt) - 1;
  int n = 0;
  for ( ; x > 0; x = x >> 1) {
    ++n;
  }
  return n-1;
}

int count_ones(char* tt) {
  int count = 0;
  for (char* p = tt; *p != '\0'; ++p) {
    if (*p == '1') {
      ++count;
    }
  }
  return count;
}

int* string_to_ints(char* tt, int const m) {
  int* T = (int*)malloc(m*sizeof(int));
  int j = 0;
  for (int i = 0; tt[i] != '\0'; ++i) {
    if (tt[i] == '1') {
      T[j] = i;
      ++j;
    }
  }
  return T;
}

void print_tt(int* tt, int m) {
  for (int i = 0; i < m; ++i) {
    printf("%d ", tt[i]);
  }
  putchar('\n');
}

plist* all_pairs_j(mpz_t* v, int l, int r, int j) {
  plist* ps = (plist*)calloc(1, sizeof(plist));
  // int m = r - l + 1; // v の長さ
  int k = 0;
  int kp = 0;

  mpz_t ZERO;
  mpz_init_set_ui(ZERO, 0);
  
  while (l+k < r) {
    while (l+k < r) {
      mpz_t res, ptr;
      mpz_init(res);
      mpz_set(res, v[l+k]);
      mpz_init_set_ui(ptr, 2);
      mpz_pow_ui(ptr, ptr, j);
      mpz_and(res, res, ptr);
      if (0 == mpz_cmp(ZERO, res)) {
	mpz_clear(ptr);
	mpz_clear(res);
	break;
      }
      ++k;
    }
    if (l+k == r) {
      return ps;
    }
    if (kp < k) {
      kp = k+1;
    }
    while (l+kp < r) {
      mpz_t res, ptr;
      mpz_init(res);
      mpz_set(res, v[l+k]);
      mpz_init_set_ui(ptr, 2);
      mpz_pow_ui(ptr, ptr, j);
      mpz_add(res, res, ptr);
      if (mpz_cmp(v[l+kp], res) >= 0) {
	mpz_clear(ptr);
	mpz_clear(res);
	break;
      }
      kp += 1;
    }
    if (r-l == kp) {
      return ps;
    }

    mpz_t res1, res2;
    mpz_init(res1);
    mpz_set(res1, v[l+k]);
    mpz_xor(res1, res1, v[l+kp]);
    mpz_init_set_ui(res2, 2);
    mpz_pow_ui(res2, res2, j+1);
    if (mpz_cmp(res1, res2) >= 0) {
      k = kp;
      mpz_clear(res2);
      mpz_clear(res1);
      continue;
    }

    mpz_t ptr;
    mpz_set(res1, v[l+kp]);
    mpz_set(res2, v[l+k]);
    mpz_init_set_ui(ptr, 2);
    mpz_pow_ui(ptr, ptr, j);
    mpz_add(res2, res2, ptr);
    if (mpz_cmp(res1, res2) == 0) {
      pair p = { k, kp };
      push_front(ps, p);
    }

    mpz_clear(ptr);
    mpz_clear(res2);
    mpz_clear(res1);

    ++k;
  }

  mpz_clear(ZERO);
  
  return ps;
}

plist* find_all_maximal_cubes(int* v, int m, int n) {
  plist* cubes = (plist*)calloc(1, sizeof(plist));

  mpz_t* S = (mpz_t*)calloc(2*m+n, sizeof(mpz_t)); // 極大でない点を管理
  mpz_t* T = (mpz_t*)calloc(2*m+n, sizeof(mpz_t)); // 極大でない点のタグを管理

  mpz_t* vv = (mpz_t*)calloc(m, sizeof(mpz_t));
  for (int j = 0; j < m; ++j) {
    mpz_init_set_ui(vv[j], v[j]);
  }

  /* for (int j = 0; j < m; ++j) { */
  /*   mpz_out_str(stdout, 2, vv[j]); */
  /*   putchar('\n'); */
  /* } */
  
  for (int j = 0; j < m; ++j) {
    plist* ps = all_pairs_j(vv, 0, m, j);

    /* plist_print_reverse(ps); putchar('\n'); */
    mpz_t ptr;
    mpz_init_set_ui(ptr, 2);
    mpz_pow_ui(ptr, ptr, j);
    for (pnode* p = ps->tail; p != NULL; p = p->prev) {
      int k = p->val.fst;
      int kp = p->val.snd;
      mpz_add(T[k], T[k], ptr);
      mpz_add(T[kp], T[kp], ptr);
    }
    mpz_clear(ptr);
    clear(ps);
  }

  int s = 0;
  int t = 0;

  mpz_t ZERO;
  mpz_init_set_ui(ZERO, 0);
  
  while (s != m) {
    if (mpz_cmp(ZERO, T[s]) == 0) {
      pair p = { 0, v[s] };
      push_front(cubes, p);
      s += 1;
    } else {
      mpz_t vs, ts;
      mpz_init(vs);
      mpz_init(ts);
      mpz_set(S[t], vv[s]);
      mpz_set(T[t], T[s]);
      t += 1;
      s += 1;
      mpz_clear(ts);
      mpz_clear(vs);
    }
  }

  mpz_t A;
  mpz_init_set_ui(A, 0);
  mpz_set_ui(S[t], 0);
  
  while (1) {
    int j = 0;
    mpz_t TT;
    mpz_init_set_ui(TT, t);
    if (mpz_cmp(S[t], TT) == 0) {
      mpz_t res;
      mpz_init(res);
      while (j < n) {
	mpz_set_ui(res, 2);
	mpz_pow_ui(res, res, j);
	mpz_and(res, res, A);
	if (mpz_cmp(ZERO, res) != 0) {
	  break;
	}
	j += 1;
      }
      mpz_clear(res);
    }
    mpz_clear(TT);
    mpz_t res;
    mpz_init(res);
    while (j < n) {
      mpz_set_ui(res, 2);
      mpz_pow_ui(res, res, j);
      mpz_and(res, res, A);
      if (mpz_cmp(ZERO, res) == 0) {
	break;
      }
      t = mpz_get_ui(S[t]) - 1;
      mpz_t ptr;
      mpz_init_set_ui(ptr, 2);
      mpz_pow_ui(ptr, ptr, j);
      mpz_sub(A, A, ptr);
      j += 1;
    }
    mpz_clear(res);
    if (j >= n) {
      break;
    } else {
      mpz_t ptr;
      mpz_init_set_ui(ptr, 2);
      mpz_pow_ui(ptr, ptr, j);
      mpz_add(A, A, ptr);
      mpz_clear(ptr);
    }

    int r = t;
    s = mpz_get_ui(S[t]);

    /* printf("%d %d: ", s, r); */
    /* for (int ii = s; ii < r; ++ii) { */
    /*   mpz_out_str(stdout, 10, S[ii]); */
    /*   putchar('\n'); */
    /* } */
    /* putchar('\n'); */
    plist* ps = all_pairs_j(S, s, r, j);
    /* plist_print_reverse(ps); putchar('\n'); */
    for (pnode* p = ps->tail; p != NULL; p = p->prev) {
      int k = p->val.fst;
      int kp = p->val.snd;

      mpz_t res1, res2;
      mpz_init_set(res1, T[s+k]);
      mpz_and(res1, res1, T[s+kp]);
      mpz_init_set_ui(res2, 2);
      mpz_pow_ui(res2, res2, j);
      mpz_sub(res1, res1, res2);
      if (mpz_cmp(ZERO, res1) == 0) {
	pair q = { mpz_get_ui(A), mpz_get_ui(S[s+k]) };
	push_front(cubes, q);
      } else {
	t += 1;
	mpz_set(S[t], S[s+k]);
	mpz_set(T[t], res1);
      }
      mpz_clear(res2);
      mpz_clear(res1);
    }
    t += 1;
    mpz_t rr;
    mpz_init_set_ui(rr, r+1);
    mpz_set(S[t], rr);
    mpz_clear(rr);
    
    clear(ps);
  }

  mpz_clear(A);
  mpz_clear(ZERO);
  
  for (int j = 0; j < m; ++j) {
    mpz_clear(vv[j]);
  }
  free(vv);
  for (int j = 0; j < 2*m+n; ++j) {
    mpz_clear(T[j]);
    mpz_clear(S[j]);
  }
  free(T);
  free(S);
  
  return cubes;
}

char* ab2str(int const a, int const b, int const n) {
  char* s = (char*)malloc((n+1)*sizeof(char));
  for (int i = 0; i < n; ++i) {
    if (0 != (a & (1 << i))) {
      s[n-i-1] = '*';
    } else if (0 != (b & (1 << i))) {
      s[n-i-1] = '1';
    } else {
      s[n-i-1] = '0';
    }
  }
  s[n] = '\0';
  return s;
}

void print_cube(pair p, int const n) {
  char* s = ab2str(p.fst, p.snd, n);
  /* printf("%s", s); */
  printf("%s (%d, %d)", s, p.fst, p.snd);
  free(s);
}

void print_cubes(plist* cubes, int const n) {
  for (pnode* c = cubes->tail; c != NULL; c = c->prev) {
    print_cube(c->val, n);
    putchar('\n');
  }
}

void print_cube2(pair p, int const n) {
  char* s = ab2str(p.fst, p.snd, n);
  /* printf("%s", s); */
  printf("%s", s);
  free(s);
}

void print_cubes2(int set, plist* cubes, int const n) {
  int i = 0;
  for (pnode* c = cubes->tail; c != NULL; c = c->prev) {
    if ((set & (1 << i)) != 0) {
      print_cube2(c->val, n);
      putchar('\n');
    }
    i += 1;
  }
}

int pop_count(int n) {
  int res = 0;
  while (n > 0) {
    if ((n & 1) == 1) {
      ++res;
    }
    n = n >> 1;
  }
  return res;
}

void print_set(pair p, int const n) {
  int a = p.fst;
  int b = p.snd;
  if (a == 0) {
    printf("%d", b);
    return ;
  }
  int const m = pop_count(a);
  int* A = (int*)malloc(m * sizeof(int));
  int i = 0;
  for (int j = 0; j < n; ++j) {
    if ((a & (1 << j)) != 0) {
      A[i] = j;
      i += 1;
    }
  }
  for (int i = 0; i < (1 << m); ++i) {
    int add = 0;
    for (int j = 0; j < m; ++j) {
      if ((i & (1 << j)) != 0) {
	add += (1 << A[j]);
      }
    }
    printf("%d ", b + add);
    /* printf("%d ", add); */
  }
  free(A);
}

void print_sets(plist* cubes, int const n) {
  for (pnode* c = cubes->tail; c != NULL; c = c->prev) {
    print_set(c->val, n);
    putchar('\n');
  }
}

int* get_set(pair p, int const n) {
  int a = p.fst;
  int b = p.snd;
  if (a == 0) {
    int* S = (int*)malloc(sizeof(int));
    S[0] = b;
    return S;
  }
  int const m = pop_count(a);
  int* A = (int*)malloc(m * sizeof(int));
  int i = 0;
  for (int j = 0; j < n; ++j) {
    if ((a & (1 << j)) != 0) {
      A[i] = j;
      i += 1;
    }
  }
  int* S = (int*)malloc(m * sizeof(int));
  for (int i = 0; i < (1 << m); ++i) {
    int add = 0;
    for (int j = 0; j < m; ++j) {
      if ((i & (1 << j)) != 0) {
	add += (1 << A[j]);
      }
    }
    S[i] = b+add;
    /* printf("%d ", b + add); */
  }
  free(A);
  return S;
}

int* converts_to_bits(int* T, int const m, plist* cubes, int const n) {
  int const SIZE = cubes->size;
  int* bits = (int*)malloc(SIZE * sizeof(int));
  int i = 0;
  for (pnode* ptr = cubes->tail; ptr != NULL; ptr = ptr->prev, ++i) {
    int* sets = get_set(ptr->val, n);
    int ones = pop_count(ptr->val.fst);
    bits[i] = 0;
    for (int j = 0; j < (1 << ones); ++j) {
      for (int k = 0; k < m; ++k) {
	if (sets[j] == T[k]) {
	  bits[i] |= (1 << k);
	}
      }
    }
    free(sets);
  }
  return bits;
}

// 整数（ビット列）で表現される集合を出力
void print_in_binary(int s, int n) {
  for (int i = 0; i < n; ++i) { 
    if (s & 1 == 1) {
      printf("1");
    } else {
      printf("0");
    }
    s = s >> 1;
  }
  putchar('\n');
}

int get_essential(int const n, int remained, int selected, int covered, int* option, int const m) {
  for (int i = 0; i < m; ++i) {
    if ((remained & (1 << i)) == 0) {
      continue;
    }
    int opt_i = option[i];
    opt_i |= covered;
    opt_i -= covered;
    for (int j = 0; j < m; ++j) {
      if ((remained & (1 << j)) == 0 || i == j) {
	continue;
      }
      opt_i |= option[j];
      opt_i -= option[j];
      if (opt_i == 0) {
	break;
      }
    }
    if (opt_i != 0) {
      return i;
    }
  }
  return -1;
}

int number_of_items(int selected, int* option, int const m) {
  int num = 0;
  for (int i = 0; i < m; ++i) {
    if ((selected & (1 << i)) != 0) {
      num += __builtin_popcount(option[i]);
    }
  }
  return num;
}

// n 個のアイテムと remained で 1 が立っている個数のオプションが存在
// remained は選択するか否かを未だ決定していないオプションの集合
// selected 選択したオプションの集合
// covered は被覆されたアイテムの集合
void mc_sub(int const n, int remained, int selected, int covered, int* option, int const m, int* res, int level) {
  /* /\********************************************\/ */
  /* for (int i = 0; i < level; ++i) printf(" "); */
  /* printf("current state: "); */
  /* print_in_binary(selected, m); */
  /* /\********************************************\/ */
  
  if ((1 << n) - 1 == covered) {
    // 解を発見したので，現在の解と比較
    // 発見した解 current のリテラル数
    int cur_nls = number_of_items(selected, option, m);
    // 現状の候補 res のリテラル数
    int res_nls = number_of_items(*res, option, m);
    if (cur_nls < res_nls) {
      *res = selected;
    }
    /* /\********************************************\/ */
    /* for (int i = 0; i < level; ++i) printf(" "); */
    /* print_in_binary(selected, m); */
    /* /\********************************************\/ */
    return;
  } else if (0 == remained) {
    // この選択 (selected) では解は見つからなかった
    return;
  }
  
  // 残っているオプションの中で必ず被覆しなければいけないもの（必須主項）がないかを探索
  while (1) {
    int essential = get_essential(n, remained, selected, covered, option, m);
    if (-1 == essential) {
      break;
    }

    /* /\********************************************\/ */
    /* for (int i = 0; i < level; ++i) printf(" "); */
    /* printf("essential = %d: ", essential); */
    /* for (int i = 0; i < level; ++i) printf(" "); */
    /* print_in_binary(option[essential], n); */
    /* /\********************************************\/ */
    
    remained -= (1 << essential);
    covered  |= option[essential];
    selected |= (1 << essential);
  }
  
  // 必須主項を選択して被覆が完了しているかを確認
  if ((1 << n) - 1 == covered) {
    // 解を発見したので，現在の解と比較
    // 発見した解 current のアイテム数
    int cur_nls = number_of_items(selected, option, m);
    // 現状の候補 res のアイテム数
    int res_nls = number_of_items(*res, option, m);
    if (cur_nls < res_nls) {
      *res = selected;
    }
    /* /\********************************************\/ */
    /* for (int i = 0; i < level; ++i) printf(" "); */
    /* print_in_binary(selected, m); */
    /* /\********************************************\/ */
    return;
  }

  // 残っているオプションの中からオプションを選択して再帰呼び出し
  for (int i = 0; i < m; ++i) {
    if ((remained & (1 << i)) == 0) {
      continue;
    }
    // オプション i を選択
    remained -= (1 << i);
    int covered_tmp = covered;
    covered  |= option[i];
    selected |= (1 << i);
    // 再帰呼出し
    mc_sub(n, remained, selected, covered, option, m, res, level+1);
    // オプション i を選択したのを undo
    remained += (1 << i);
    covered = covered_tmp;
    selected -= (1 << i);
  }
}


// 最小被覆
// n 個のアイテムと m 個のオプション (
int mc(const int n, int* option, const int m) {
  /* for (int i = 0; i < m; ++i) { */
  /*   print_in_binary(option[i], n); */
  /* } */
  
  int remained = (1 << m) - 1;
  int selected = 0;
  int covered = 0;
  int res = (1 << m) - 1;
  mc_sub(n, remained, selected, covered, option, m, &res, 0);
  
  return res;
}

int main()
{
  char tt[MAX];
  fgets(tt, MAX, stdin);

  int const n = get_number_of_variables(tt);  
  int const m = count_ones(tt);

  int* T = string_to_ints(tt, m);
  /* print_tt(T, m); */

  /* printf("m = %d\n", m); */
  
  plist* cubes = find_all_maximal_cubes(T, m, n);

  /* plist_print_reverse(cubes); putchar('\n'); */
  /* /\* */
  /*  * (a, b) */
  /*  * a は * の箇所 */
  /*  * b は 1 が立っている箇所 */
  /*  *\/ */
  
  /* print_cubes(cubes, n); */
  /* print_sets(cubes, n); */
  /* printf("%d %d\n", m, cubes->size); */
  /* printf("==========\n"); */

  int* bits = converts_to_bits(T, m, cubes, n);
  int res = mc(m, bits, cubes->size);

  /* print_in_binary(res, cubes->size); */
  
  print_cubes2(res, cubes, n);
  
  free(bits);
  clear(cubes);
  free(T);
  
  return 0;
}