Sunday, September 21, 2008

Number Split

topCoder: NumberSplit

By Editorial:

The first thing to notice here is that the numbers in every step become smaller. Let's take an example and say we have a five digit number of the form abcde (a, b, c, d and e represent the decimal digits of the number), which we split to produce the successor: ab * c * de. Here it is always c < 10 and de < 100, so for the successor we have: ab * c * de < ab * 10 * 100 = ab000 ≤ abcde. Similar with any other numbers and splittings.

Now we need a way to compute all possible successors, given a given number. For this we can use the following recursive pseudo-code:

generateSuccessors(int multiplier, int n) {
add (multiplier * n) to the set of successors
for (i = 10; i <= n; i *= 10) {
generateSuccessors(multiplier * (n / i), n % i);
}
}


We initialize the set of successors to an empty set, and call generateSuccessors(1, n) (where n is the number, for which we want to find the successors). Finally, we remove n from the generated set (since this is not a successor of n itself, we need to split the given number to at least two parts for a successor to be valid).


To compute the longest possible sequence, starting with the given start, generate all successors of start as described above, and for each n in the successor set compute recursively longestSequence(n). The return value is the maximum of all computed values + 1 (if start is a single digit number, the successors set will be empty, and we return 1). Note that this would not work if loops in the sequence were possible, but since each successor is smaller then the original number, this can not happen. In order to avoid a timeout, we need to memorize in a buffer all values already computed.


Alternatively, we can use dynamic programming, by initializing longest[i] = 1 for all single digit numbers i, and computing longest[i] from i = 10 up to i = start by adding 1 to the maximum of longest[j] for all j in the successor set of i. This works, since all successors of i are smaller than i. With this solution, we compute the longest sequence for more numbers than actually needed (many numbers between 1 and start can not be reached from start) but with the low constraints this solution was within the time limit.


my dp code with memoization (note: topCoder does not have itoa)


void itoa(int n, char* buff) {
sprintf(buff, "%d", n);
}
// memoization array
int a[999999];
int solve(string s) {
int ss = atoi(s.c_str());
if(a[ss] != -1) return a[ss];
char buff[10];
if(s.length() == 1) return 1;
int res = 0;
for(int i = 1; i < s.length(); ++i) {
// single split
int left = atoi(s.substr(0,i).c_str());
int right = atoi(s.substr(i).c_str());
itoa(left*right, buff);
checkmax(res, 1+solve(string(buff)));
string rights = s.substr(i);
if(rights.length() <= 1) continue;

// this part does double split
for(int j = 1; j < rights.length(); ++j) {
int lleft = atoi(rights.substr(0,j).c_str());
int lright = atoi(rights.substr(j).c_str());
itoa(left * lleft * lright, buff);
checkmax(res, 1+solve(string(buff)));
}
}
return a[ss] = res;
}

int NumberSplit::longestSequence(int start) {
char buff[10];
itoa(start,buff);
memset(a, -1, sizeof(a));
string s(buff);
return solve(buff);
}

Floyd-Warshall: All-Pairs Shortest Path

By wiki, Floyd-Warshall is a graph analysis algorithm for finding shortest paths in a weighted, directed graph. A single execution of the algorithm will find the shortest paths between all pairs of vertices. The Floyd–Warshall algorithm is an example of dynamic programming.

 1 /* Assume a function edgeCost(i,j) which returns the cost of the edge from i to j
2 (infinity if there is none).
3 Also assume that n is the number of vertices and edgeCost(i,i)=0
4 */
5
6 int path[][];
7 /* A 2-dimensional matrix. At each step in the algorithm, path[i][j] is the shortest path
8 from i to j using intermediate values in (1..k-1). Each path[i][j] is initialized to
9 edgeCost(i,j).
10 */
11
12 procedure FloydWarshall ()
13 for k: = 0 to n − 1
14 for each (i,j) in (0..n − 1)
15 path[i][j] = min ( path[i][j], path[i][k]+path[k][j] );

For numerically meaningful output, Floyd-Warshall assumes that there are no negative cycles (in fact, between any pair of vertices which form part of a negative cycle, the shortest path is not well-defined because the path can be infinitely small). Nevertheless, if there are negative cycles, Floyd–Warshall can be used to detect them. A negative cycle can be detected if the path matrix contains a negative number along the diagonal. If path[i][i] is negative for some vertex i, then this vertex belongs to at least one negative cycle.

Applications:


  • Shortest paths in directed graphs (Floyd's algorithm).
  • Transitive closure of directed graphs (Warshall's algorithm). In Warshall's original formulation of the algorithm, the graph is unweighted and represented by a Boolean adjacency matrix. Then the addition operation is replaced by logical conjunction (AND) and the minimum operation by logical disjunction (OR).
  • Optimal routing. In this application one is interested in finding the path with the maximum flow between two vertices. This means that, rather than taking minima as in the pseudocode above, one instead takes maxima. The edge weights represent fixed constraints on flow. Path weights represent bottlenecks; so the addition operation above is replaced by the minimum operation.
  • Testing whether an undirected graph is bipartite
  • Minmax / Maxmin Distance

    Shortest Path Problems:

    int Floyd_Warshall (int n) {
    for(int k = 0; k < n; ++k)
    for(int i = 0; i < n; ++i)
    for(int j = 0; j < n; ++ j)
    G[i][j] = min(G[i][j], G[i][k] + G[k][j]);
    }


    Transitive Closure/Hull

    for(int k = 0; k < N; ++k)
    for(int i = 0; i < N; ++i)
    for(int j = 0; j < N; ++j)
    G[i][j] = G[i][j] || (G[i][k] && G[k][j]);


    Minmax/ Maxmin distance

    void minmax(int n) {


    for(int k = 0; k < n; ++k)
    for(int i = 0; i < n; ++i)
    for(int j = 0; j < n; ++j)
    G[i][j] = min(G[i][j], max(G[i][k], G[k][j]));
    }



      void maxmin(int n) {
      for(int k = 0; k < n; ++k)
      for(int i = 0; i < n; ++i)
      for(int j = 0; j < n; ++j)
      G[i][j] = max(G[i][j], min(G[i][k], G[k][j]));
      }


    • UVa 534: Frogger

    • UVa 10048: Audiophobia

    • UVa 544: Heavy Cargo

    • UVa 10099: The Tourist Guide

  • Subset Sum Problem

    The problem is this: given a set of integers, does the sum of some non-empty subset equal exactly zero?The problem is NP-Complete. wiki

    An equivalent problem is this: given a set of integers and an integer s, does any non-empty subset sum to s? Subset sum can also be thought of as a special case of the knapsack problem. One interesting special case of subset sum is the partition problem, in which s is half of the sum of all elements in the set.

    topCoder: PayBill Editorial

    There are two basic approaches to this problem. The first, which inevitably fails due to timeout on larger test cases, is to try obtain the sum by either including or excluding the first element, and then calling itself recursively with the remainder of the set. Unfortunately, with the maximum 50 people, this is over 1 quadrillion operations to perform.

    Instead, a more clever, dynamic programming approach needs to be used. Since each item can only be up to 10,000, we know that the sum cannot be more than 500,000. So, we simply need a boolean array of 500,000 elements, where the i-th element represents whether or not we can reach the total i with some subset of the original values. In code it looks something like this:

    boolean[] canTotal = new boolean[500001];
    canTotal[0] = true;
    for (int i = 0; i < meals.length; i++)
    for (int j = totalMoney; j >= meals[i]; j--)
    if (canTotal[j - meals[i]]) canTotal[j] = true;

    topCoder: LoadBalancing Editorial


    int minTime(vector<int> chunkSizes)
    {
    vector<int> dp(204801, 0);
    dp[0] = 1;
    int total = 0;
    for (int i=0; i<(int)chunkSizes.size(); i++)
    {
    total += chunkSizes[i]/1024;
    for (int j=204800; j>=0; j--)
    if (dp[j] == 1)
    dp[j+chunkSizes[i]/1024] = 1;
    }
    for (int i=(total+1)/2; true; i++)
    if (dp[i] == 1)
    return 1024 * i;
    }

    UVa 562: Dividing Coins



    跟上面一题如出一辙,就是平分问题。这里用了bitvector

    const int MAX = 50000;
    const unsigned int N = 1563;
    unsigned int a[N+1];

    int coin[101];
    int main() {
    int ndata;
    scanf("%d ", &ndata);
    while(ndata--) {
    int n;
    scanf("%d ", &n);

    for(int i = 0; i < n; ++i) {
    scanf("%d ", &coin[i]);
    }

    int total = std::accumulate(coin, coin+n, 0);

    memset(a, 0, sizeof(int)*(total>>SHIFT));
    SET(a, 0);
    for(int i = 0; i < n; ++i) {
    for(int j = total; j >= 0; --j)
    if(TEST(a, j))
    SET(a, j+coin[i]);
    }

    for(int i = (total+1)/2; true; ++i) {
    if(TEST(a, i)) {
    printf("%d\n", abs(2*i-total));
    break;
    }
    }
    }
    }

    UVa 11517 Exact Change



    这道要多记一个coin的数目,不能用bitvector了

            a[0] = 1;
    for(int i = 0; i < n; ++i) {
    for(int j = max-1; j >= 0; --j)
    if(a[j]) {
    int t = j+ bill[i];
    if(t > max) continue;
    if(a[t] != 0)
    a[t] = a[t] < a[j]+1 ? a[t]: a[j]+1;
    else
    a[t] = a[j]+1;
    }
    }
    for(int i = price; i <= MAX; ++i) {
    if(a[i]) {
    printf("%d %d\n", i, a[i]-1);
    break;
    }
    }

    Saturday, September 20, 2008

    memo[i][j] = min ( memo[i][k] + memo[k+1][j] )

    用类似这一recurrence的题非常多,在这汇总一下

      topCoder: QuickSums

      INPUT: i, j, sum
      for all i <= k <= j take out substring(i,k) from sum try d[k+1][j] = solve(substring(k+1,j), sum - substring(i,k)); end for d[i][j] = min (1+d[k+1][j]) k=i...j

      边缘情况注意一下,如果d[i][j]整体已经等于sum,直接返回,不用1+res

      #define INF 999999999
      int memo[12][12][102];
      void checkmin(int &a, int b) {
      if(b < a) a = b;
      }
      int solve(string numbers, int i, int j, int s) {
      if(memo[i][j][s] >= 0) return memo[i][j][s];
      if(i == j) {
      if (numbers[i]-'0' == s) return memo[i][j][s] = 0;
      else return memo[i][j][s] = INF;
      }
      // the whole word match!
      if(atoi(numbers.substr(i, j-i+1).c_str()) == s) return memo[i][j][s] = 0;
      int res = INF;
      for(int k = i; k < j; ++k) {
      int left = atoi(numbers.substr(i,k-i+1).c_str());
      if(left <= s){
      //cout << numbers.substr(k+1, j-k+1) << "->" << s-left << endl;
      int sol = 1+solve(numbers, k+1, j, s-left);
      //cout << numbers.substr(k+1, j-k+1) << " (" << sol <<")" << endl;
      checkmin(res, sol);
      }
      }
      return memo[i][j][s] = res;
      }

      int QuickSums::minSums(string numbers, int sum) {
      memset(memo, -1, sizeof(memo));
      int res = solve(numbers, 0, numbers.length()-1, sum);
      if(res == INF) return -1;
      return res;
      }


      topCoder: ShortPalindromes


      “ shortest(base)
      if base is already a palindrome then
      return base
      if base has the form A...A then
      return A + shortest(...) + A
      if base has the form A...B then
      return min(A + shortest(...B) + A,
      B + shortest(A...) + B)
      ”
      string memo[26][26];
      string solve(string base, int i, int j) {
      //cout << base.substr(i,j-i+1) << endl;
      if(i>j) return string("");
      if(i == j) return memo[i][j]=base[i];
      if(memo[i][j] != "") return memo[i][j];

      string res;
      if(base[i] == base[j]) {
      res = base[i] + solve(base, i+1, j-1)+base[i];
      return memo[i][j]=res;
      }

      string pre = base[j] + solve(base, i, j-1) + base[j];
      string suf = base[i] + solve(base, i+1, j) + base[i];

      if(pre.length() != suf.length()) {
      if(pre.length() < suf.length()) res = pre;
      else res = suf;
      }
      else {
      if(pre < suf) res = pre;
      else res = suf;
      }
      return memo[i][j]=res;
      }

      string ShortPalindromes::shortest(string base) {

      int L = base.length();
      for(int i = 0; i < L; ++i)
      for(int j = 0; j < L; ++j)
      memo[i][j]="";
      return solve(base, 0, base.length()-1);
      }

      topCoder: TreePlanting


      long long memo[62][62][62];
      int N;
      long long solve(long long status, int i, int j, int f) {

      if(memo[i][j][f] != -1LL) return memo[i][j][f];

      if(f == 0 || (i == j && f == 1)) return memo[i][j][f] = 1LL;
      if(i > j) return memo[i][j][f] = 0LL;
      if(i == j && f > 1) return memo[i][j][f] = 0LL;

      long long res = 0;
      for(int k = i; k <= j; ++k) if(status & 1LL<<k) {

      res += solve(status ^ (1LL<<k), k+2, j, f-1);
      }
      return memo[i][j][f] = res;
      }

      long long TreePlanting::countArrangements(int total, int fancy) {
      memset(memo, -1, sizeof(memo));
      N = total;
      long long status = (1LL<<total)-1;
      return solve(status, 0, total-1, fancy);
      }
      Editorial solution:

      Here, we need to get a little creative to come up with a workable solution. Consider placing f fancy trees along a total of n locations, so that no two are adjacent. Call our function, C, the number of ways to do this. At our first location, we can either plant a fancy tree, or not plant a fancy tree. If we plant a fancy tree, then we can't plant a fancy tree in the second spot, and hence will have n-2 locations in which to plant the remaining f-1 fancy trees. If we don't plant a fancy tree in the first spot, then we have n-1 locations in which to place all f fancy trees. This gives us our recursive formula, C(n, f) = C(n - 2, f - 1) + C(n - 1, f). Our starting values are of course, C(1, 1) = 1, and C(0, 0) = 1.

      However, with the problem constraints as they are, a simple recursive function by itself is not sufficient. We either need to implement memoization, whereby we store the values of previous recursive calls, to avoid recalculating them repeatedly, or we can use dynamic programming, as shown here:

      public long countArrangements(int total, int fancy) {
      long[][] count = new long[total + 1][fancy + 1];
      count[0][0] = 1;
      count[1][1] = 1;
      for (int i = 0; i < = total; i++)
      for (int j = 0; j < = fancy; j++) {
      if (i > 0)
      count[i][j] += count[i - 1][j];
      if (i > 1 && j > 0)
      count[i][j] += count[i - 2][j - 1];
      };
      return count[total][fancy];
      }

      topCoder: SentenceDecomposition

      int diff(string& a, string& b) {
      int cost = 0;
      for(int i = 0;i < a.length(); ++i)
      if(a[i]!=b[i]) cost ++;
      return cost;
      }
      bool match(string& aa, string &bb) {
      sort(ALL(aa));
      return aa==bb;
      }
      int memo[52];
      //dict is a copy of validWords
      //sort_dict is a copy of dict, but each row is sorted
      vector<string> dict, sort_dict;
      #define INF 999999999
      int solve(int i, string &s) {
      //cout << s.substr(i) << endl;
      if(memo[i] != -1) return memo[i];

      if(i == s.length()) return memo[i] = 0;
      if(i > s.length()) return memo[i]=INF;

      bool found = false;
      int mincost = INF;
      for(int k = 0; k < dict.size(); ++k) {
      int len = dict[k].length();
      if(!match(s.substr(i, len), sort_dict[k])) continue;
      int ret = 0;
      ret = solve(i+len, s);
      if(ret == INF) continue;
      else {
      found = true;
      ret += diff(s.substr(i, len), dict[k]);
      }
      mincost = min(mincost, ret);
      }

      if(found) return memo[i] = mincost;
      return memo[i] = INF;
      }

      int SentenceDecomposition::decompose(string sentence, vector <string> validWords) {
      sort_dict = dict = validWords;
      for(int i = 0; i < sort_dict.size(); ++i)
      sort(ALL(sort_dict[i]));

      memset(memo, -1, sizeof(memo));

      int res = solve(0, sentence);
      if(res == INF) return -1;
      return res;
      }

      Friday, September 19, 2008

      Longest Increasing/Decreasing Subsequence

      A lot of related problems:


      two basic dp algorithms

      e.g.,
      sequence: 1,6,2,3,5,4,7
      algorithm 1: O(N^2)
      a 1 2 3 4 5 6 7
      ---------------------------
      1 6 2 3 5 4 7
      ---------------------------
      0| 1 1 1 1 1 1 1
      1| (1) 2 2 2 2 2 2
      2| 1 (2) 2 2 2 2 2
      3| 1 2 (2) 3 3 3 3
      4| 1 2 2 (3) 4 4 4
      5| 1 2 2 3 (4) 4 5
      6| 1 2 2 3 4 (4) 5 <-

      algorithm 2: O(NlogN)
      a 1 2 3 4 5 6 7
      ---------------------------
      1 6 2 3 5 4 7
      ---------------------------
      0| 1
      1| 1 6
      2| 1 2
      3| 1 2 3
      4| 1 2 3 5
      5| 1 2 3 4
      6| 1 2 3 4 7 <-
      topCoder: thePriceIsRight 
      // longest increasing sequence (O(N^2)) 
      vector<int> lis (vector<int> prices) {
      vector<int> len(prices.size(), 1); //length of the longest inc sequence up to i
      vector<int> ways(prices.size(), 1); //how many ways of the lcs up to i

      for(int i = 0; i < prices.size()-1; ++i) {
      for(int j = i+1; j < prices.size(); ++j)
      if(prices[j] > prices[i]) {
      if(len[j] < len[i] + 1) {
      len[j] = len[i] + 1;
      ways[j] = ways[i];
      }
      else if(len[j] == len[i]+1)
      ways[j]+= ways[i];
      }
      }
      int maxv = *max_element(len.begin(), len.end());
      int nways = 0;
      for(int i = 0; i < len.size(); ++i)
      if(len[i] == maxv) nways += ways[i];
      vector<int> res;
      res.push_back(maxv), res.push_back(nways);
      return res;
      }

      topCoder: Books
      // binary search 
      template<typename T1, typename T2> int bsearch(T1 &A, T2 target) {
      int u = 0;
      int v = A.size()-1;
      while(u < v) {
      int c = (u + v) / 2;
      if (target > A[c]) u=c+1;
      else v=c;
      }
      return u;
      }

      // longest non-decreasing sequence
      template<typename T> int lis (vector<T> &seq) {
      if(seq.empty()) return 0;
      vector<T> A;
      int n = seq.size();
      A.push_back(seq[0]);

      int u, v;
      for(int i = 1; i < n; ++i) {
      if(seq[i] >= A.back()) {
      A.push_back(seq[i]);
      continue;
      }
      int u = bsearch(A, seq[i]);
      //important !!! when compute non-increasing or non-decreasing seq
      if(A[u] == seq[i]) A[u+1] = seq[i];
      else A[u] = seq[i];
      }
      return A.size();
      }

      int Books::sortMoves(vector <string> titles) {
      return books.size()-lis(titles);
      }

      WordParts

      TopCoder: WordParts –> Editorial

      • prefix, suffix
      • memoization

      d[i][j] = min(d[i][k] + d[k+1][j])
      k --- substring(i,k) exists in the prefix/suffix array
      set<string> dict;
      #define INF 999999999
      int memo[51];

      void checkmin(int &a, int b) {
      if (b < a) a = b;
      }
      int solve(string &compound, int cpos) {
      //cout << compound.substr(cpos) << " " << cpos << endl;
      if(memo[cpos] >= 0) return memo[cpos];
      if(cpos == compound.length()) return memo[cpos] = 0;

      int res = INF;
      int maxlen = compound.length() - cpos;

      for(int i = 1; i <= maxlen; ++i) {
      string s = compound.substr(cpos, i);
      if(dict.count(s)) {
      checkmin(res, 1+solve(compound, cpos+i));
      }
      }
      return memo[cpos] = res;
      }

      int WordParts::partCount(string original, string compound) {
      dict.clear(); // don’t forget to reset dict
      // add the whole word
      dict.insert(original);
      // add all prefix
      for(int i = 1; i <= original.length()-1; ++i)
      dict.insert(original.substr(0,i));
      // add all suffix
      for(int i = 1; i <= original.length()-1; ++i)
      dict.insert(original.substr(i));

      memset(memo, -1, sizeof(memo));
      int res = solve(compound, 0);
      if(res == INF) return -1;
      return res;
      }

      dancing couples

      topCoder: DancingCouples editorial
      "Each subproblem is clearly defined by three variables: the set of unmatched boys, the set of unmatched girls, and the number of couples we need to make."
      d[i][j][k]  ----- the total ways to make k couples 
      | | |
      | | i couples to make
      | available girl status (bitmask)
      the i-th boy to match

      d[i][j][k] = the ways to make k couples in which we do not use the i-th boy +
      the ways to make k couples in which we use the i-th boy
      = d[i-1][j][k] + sum (d[i-1][j^(1 << k)][k-1])
      {k|A[k][j]='Y'}
      对比上一贴:

      cnt [v][s] = the ways in which z6 = 0 + the ways in which z6 > 0
      = cnt[v-1][s] + cnt[v][s-v]
      The code:
      int dp[12][1025][12];
      int solve(int nboy, int available_girl, int K) {
      if(nboy < K || count_bit(available_girl) < K) return 0;
      if(K == 0) return 1;
      if(dp[nboy][available_girl][K] >= 0) return dp[nboy][available_girl][K];

      int res = solve(nboy-1, available_girl, K);

      for(int i = 0; i < cd[0].size(); ++i) {
      if(cd[nboy-1][i] == 'Y' && available_girl & (1<<i))
      res += solve(nboy-1, available_girl ^ (1<<i), K-1);
      }
      return dp[nboy][available_girl][K] = res;
      }

      knapsack problem

      0/1 Knapsack Problem

      n items (U = {u1, u2, … , u_n} ) need to be packed in a knapsack of size C. Each item has value v_i and size s_i.We want to find a subset of U such that

      sum_ vi is maximized subject to the constraint

      sum_si <= C

      let V[i][j] denote the value obtained by filling a knapsack of size j with items taken from the first i items in an optimal way.
      • V[i][j] = 0 (if i = 0 or j = 0)
      • V[i][j] = V[i-1][j] ( if j < si )
      • V[i][j] = max ( V[i-1][j], V[i-1][j-si] + vi ) ( if i > 0 and j >= si )

      相关的dp题
      SuperSale:: 标准knapsack, 每人来一轮knapsack
      #include <iostream>
      #include <vector>
      #include <algorithm>

      using namespace std;

      int C[1001][31];

      int Vi[1001], Wi[1001];
      int MW[31];

      void knapsack(int N, int MaxW) {
          for(int i = 0; i < N; ++i)
              C[i][0] = 0;
          for(int w = 0; w <= MaxW; ++w)
              C[0][w] = 0;
          for(int i = 1; i <= N; ++i) {
              for( int w = 1; w <= MaxW; ++w) {
                  if(Wi[i] > w)
                      C[i][w] = C[i-1][w];
                  else
                      C[i][w] = max(C[i-1][w], C[i-1][w-Wi[i]] + Vi[i]);
              }
          }
      }

      int main() {

          int T; cin >> T;

          for(int t = 0; t < T; ++t) {
              int N; cin >> N;
              for(int i = 1; i <= N; ++i) {
                  cin >> Vi[i] >> Wi[i];
              }

              int G; cin >> G;
              for(int i = 0; i < G; ++i)
                  cin >> MW[i];
              int total = 0;
              int maxw = *max_element(MW, MW + G);
              knapsack(N, maxw);
              for(int i = 0; i < G; ++i)
                  total += C[N][MW[i]];
              cout << total << endl;
          }

          return 0;
      }


      topcoder editorial 对第一题有很详细的讲解,中心意思是简化问题先,将原问题转化为一个更易表达和编程的形式。最后是这个形式:
      0 ≤ z1, ..., z6
      z1 + 2z2 + ... + 6z6 ≤ 6M-21

      "We can now count the number of valid sequences as follows. There are two types of valid sequences: Those where z6>0 and those where z6=0. In the first case, we can subtract 1 from z6, and get the same problem with the right side of the inequality smaller by 6. In the second case, we get a similar problem with only 5 variables.

      This can be formulated as a recurrence relation. Let cnt[v][s] be the number of ways in which we can set variables z1 to zv so that the total (z1 + ... + vzv) does not exceed s. Then we have: cnt[v][s] = cnt[v-1][s] + cnt[v][s-v]. "

      第二个题即:
      0≤z1≤a1,..., 0≤z6≤a6
      z1 + 2z2 + ... + 6z6 = (a1 + a2 + ... + a6)/2
      如果没有a1,...,a6的限制,跟题1几乎一样了,只要把"+"改成"||"就行了。
      cnt[v][s] = cnt[v-1][s] || cnt[v][s-v].

      但跟题1的一点不同是z1,...,z6有上限,所以dp过程中我们记录了zv已经用了几次。而且这题只要知道能不能"="(平分)即可,我们用0表示否,1表示可以(zv一次没用),2(zv用了一次),..., i(zv用了i-1次了),...
              int maxsum = 0;
      for(int v = 1; v <= 6; ++v) {
      maxsum += v * a[v];
      if (maxsum > sum) maxsum = sum;
      cnt[v][0] = 1; // zero can be partitioned
      for(int s = 0; s <= maxsum; ++s) {
      cnt[v][s] = cnt[v-1][s] > 0;
      if(cnt[v][s] == 1) continue;
      if(s >= v && cnt[v][s-v]) {
      if(cnt[v][s-v] - 1 < a[v]) cnt[v][s] = cnt[v][s-v]+1;
      }
      }
      }

      下面完整的code也用了题1的技巧,用v%2省了2/3的memory,速度也由63ms降为1ms以下。因为dp每次只用上一轮的结果,e.g.,算cnt[5][s]的时候只用cnt[4][s]。
      #include <iostream>
      #define N 20000

      int cnt[2][6*N/2+1];

      int main() {
      int a[7] = {0};
      int t = 1;
      while(scanf("%d %d %d %d %d %d",
      &a[1],&a[2],&a[3],&a[4],&a[5],&a[6])!= EOF &&
      (a[1] || a[2] || a[3] || a[4] || a[5] || a[6])) {

      printf("Collection #%d:\n", t++);
      int sum = 0;
      for(int i = 1; i <= 6; ++i)
      sum += i * a[i];

      // odd sum
      if(sum & 1) {
      printf("Can't be divided.\n");
      printf("\n");
      continue;
      }
      sum = sum >> 1; // divided by 2
      memset(cnt, 0, sizeof(cnt));

      cnt[0][0] = 1;

      int maxsum = 0;
      for(int v = 1; v <= 6; ++v) {
      maxsum += v * a[v];
      if (maxsum > sum) maxsum = sum;
      cnt[0][0] = 1; // zero can be partitioned
      for(int s = 0; s <= maxsum; ++s) {
      cnt[v%2][s] = cnt[1-v%2][s] > 0;
      if(cnt[v%2][s] == 1) continue;
      if(s >= v && cnt[v%2][s-v]) {
      if(cnt[v%2][s-v] - 1 < a[v]) cnt[v%2][s] = cnt[v%2][s-v]+1;
      }
      }
      }

      if(cnt[0][sum]) printf("Can be divided.\n");
      else printf("Can't be divided.\n");
      printf("\n");
      }
      return 0;
      }
      来一个例子:

      -----------------------
      z1 z2 z3 z4 z5 z6
      -----------------------
      1 1 2 1 1 1
      -----------------------

      sum = 24, target = sum/2 = 12
      -------------------------------------
      s |[0] [1] [2] [3] [4] [5] [6] [7] [8] [9] [10][11][12][13]
      ---------+----------------------------------------------------------
      cnt[1][s]| 1 2 0
      --------------------------------------------------------------------
      cnt[2][s]| 1 1 2 2 0
      --------------------------------------------------------------------
      cnt[3][2]| 1 1 1 1 2 2 2 3 3 3 0
      --------------------------------------------------------------------
      cnt[4][2]| 1 1 1 1 1 1 1 1 1 1 2 2 2 0
      --------------------------------------------------------------------
      cnt[5][2]| ......
      --------------------------------------------------------------------


      又碰到一道类似的,也放在这了
      topCoder: windowWasher
      0≤z1,z2,...,zn
      z1 + z2 + ... + zn = width of the wall (每种worker的人数之和等于wall width)
      minimize max(T1 x z1, T2 x z2, ..., Tn x Zn) (最小化总时间,大家可以并行干活)
      dp recurrence:
      cnt[v][s] = min(cnt[v-1][s], max(cnt[v-1][s-n], T[v]*n))
      n=1..s
      = min(max(cnt[v-1][s-n], T[v]*n))
      n=0..s
      启用第v个worker时,两种情况:
      1. 不用这个worker,只用前v-1个worker
      2. 用这个worker 1次,2次, ... s次
      只比前两题多一层循环,同样可以用% 2的方法,这里没用了。
      int cnt[51][1001];

      int WindowWasher::fastest(int width, int height, vector <int> washTimes) {
        int N = washTimes.size();
        vector<int> T(N+1);
        for(int i = 1; i <= N; ++i)
          T[i] = height * washTimes[i-1];
        memset(cnt, -1, sizeof(cnt));

        // use only the first worker
        for(int s = 0; s <= width; ++s)
          cnt[1][s] = s*T[1];

        // start from only one worker, add worker one by one
        for(int v = 2; v <= N; ++v) {
          cnt[v][0] = 0; // no work! no time needed!
          for(int s = 1; s <= width; ++s) {
            cnt[v][s] = cnt[v-1][s]; // do not use worker v, use v-1 worker only
            // worker v works on 1,2,...,s columns
            for(int n = 1; n <= s; ++n)
              checkmin(cnt[v][s], max(cnt[v-1][s-n], T[v]*n));
          }
        }
        return cnt[N][width];
      }

      Thursday, September 18, 2008

      Matrix Chain Multiplication

      经典教科书题
      N个矩阵,N+1个维数

      Matrix: 0 ...... N-1
      -----------+---------------------------------------------------------
      Matrix dim:| (M0 x M1) x (M1 x M2) x (M2 x M3) x ... x (M_N-1 x M_N)
      -----------+---------------------------------------------------------
      Matrix | 0 1 2 N-1
      ---------------------------------------------------------------------

      The i-th matrix's dimension: M_(i-1) x M_i
      第i到j个矩阵相乘的代价

      d[i][i] = 0
      d[i][j] = min(d[i][j], d[i][k]+d[k+1][j]+ M_(i-1)*M_(k)*M_(j))
      int C[N+1][N+1];
      int K[N+1][N+1];
      void mcm(vector<int> &M, int n) {
      for(int i = 1; i <= n; ++i)
      C[i][i] = 0;
      for(int d = 2; d <= n; ++d) {
      for(int i = 1; i <= n-d+1; ++i) {
      int j = i + d - 1;
      C[i][j] = INF;
      for(int k = i; k <= j-1; ++k) {
      int tmp = C[i][k] + C[k+1][j] + M[i-1] * M[k] * M[j];
      if( tmp < C[i][j]) {
      K[i][j] = k;
      C[i][j] = tmp;
      }
      }
      }
      }
      }
      C[1][N]中存着从1到N连乘最小的代价,从K[][]可以得到路径。
      void print_opt(int i, int j) {
      if(i == j)
      cout << "A" << i;
      else {
      cout << "(";
      print_opt(i, K[i][j]);
      cout << " x ";
      print_opt(K[i][j]+1, j);
      cout << ")";
      }
      }
      调用时传入mcm()一个长度为N+1的数组,和N(不是N+1!)
      vector<int> M(n+1);
      mcm(M, n);
      下面是memoization的版本,编程简单
      int lookupMcm(vector<int> &M, int i, int j) {
      if(C[i][j] != INF) return C[i][j];
      if(i == j)
      return C[i][j] = 0;
      for(int k = i; k < j; ++k) {
      int q = lookupMcm(M, i, k) + lookupMcm(M, k+1, j) + M[i-1] * M[k] * M[j];
      if( C[i][j] > q)
      C[i][j] = q;
      }
      return C[i][j];
      }

      调用前将C[][]初始化为INF
      for(int i = 1; i <= n; ++i)
      for(int j = 1; j <= n; ++j)
      C[i][j] = INF;

      printf("%d\n", lookupMcm(M, 1, n));
      Similarly, this technique can be used to solove the following problem: (in today's interview (oct 21) )

      find the length of the longest regular brackets sequence that is a subsequence of s

      主要思想就是从中间出发,往两边发展,

      如果已和dp[i][j]则考虑 dp[i-1][j+1]

      if  s[i-1] match s[j+1], then dp[i-1][j+1] = dp[i][j]

      else  dp[i-1][j+1] = max (dp[i-1][k] + dp[k+1][j+1])

      #include <cstdio>
      #include <cstring>

      const int MAX = 101;

      // dynamic programming 2d array
      int dp[MAX][MAX];

      // return true if find matching brackets
      inline bool match(char a, char b) {
      if(a == '(' && b == ')') return true;
      if(a == '[' && b == ']') return true;
      return false;
      }

      int main() {
        char buff[MAX];

      // terminate until "end"
      while (gets(buff) && buff[0] != 'e') {
      int n = strlen(buff);
      memset(dp, 0, sizeof(dp));
      for (int i = 0; i < n; ++ i)
      if (match(buff[i],buff[i+1]))
      dp[i][i+1] = 2;

      for (int k = 2; k < n; ++ k) {
      for (int i = 0; i < n; ++ i) {
      if (k + i < n) {
      if (match(buff[i], buff[i+k]))
      dp[i][i+k] = dp[i+1][i+k-1] + 2;
      for (int j = i; j < i + k; ++ j) {
      if (dp[i][j] + dp[j+1][i+k] > dp[i][i+k])
      dp[i][i+k] = dp[j+1][i+k] + dp[i][j];
      }
      }
      }
      }
      printf("%d\n", dp[0][n-1]);
      }
      return 0;
      }

      Wednesday, September 17, 2008

      Triangle War II (ZJU Sep 2008 contest)

      ZOJ Problem Set - 3038

      一开始以为是lights out这种puzzle,仔细看发现只能按unpushed button, 这是与lights out 类问题严重的不同,不是finite field。只好上brute force, 最多状态,2^21还是可行的,先试了DFS,发现stack overflow, 然后用BFS搞定,发生了一点小问题,visited忘了reset...

      这道题里用bit operation很方便的。
      • 首先把adj matrix表示成N个bit vector,每个bit vector是从一个结点出发的所有边。
      • 同样将当前图的状态也表示成一个N位的bit vector (status, 1 表示unpushed, 0表示pushed)。
      • 用bit vector i 同当前graph结点的状态 status 异或就直接得到push button i后的结点状态。
      下面的code是查看当前graph结点i 是否是unpushed, 如果是则按下它得到之后的状态。
      if(status &(1 << i))       
      status = status ^ adj[i];
      #include <queue>
      #include <ctime>
      using namespace std;

      // triangle size
      const int MAX_LEVEL = 6;

      // node # in the triangle
      const int MAX_N = (MAX_LEVEL+1)*(MAX_LEVEL)/2;

      short visited[(1<<MAX_N)+1];
      // adjancency matrix
      int adj[MAX_N];
      int n; // # of levels
      int maxon; // max # of unpushed node

      void construct_adj(int N) {
      memset(adj, 0, sizeof(adj));
      for(int level = 1; level <= n; ++level) {
      int i = level*(level-1)/2; // start point of this level
      for(int j = 0; j < level; ++j) {
      int a = i+j;
      adj[a] |= (1<<a); // self connected
      if(level < n) {
      adj[a] |= 1<<(a+level); // connect with left child
      adj[a] |= 1<<(a+level+1); // connect with right child
      }
      if(j < level-1) adj[a] |= 1<<(a+1); // connect with right sibling
      if(j > 0) adj[a] |= 1<<(a-1); // connect with left sibling
      }
      }
      // sysmetric filling
      for(int i = 1; i < N; ++i)
      for(int j = 0; j < i; ++j) {
      int t = adj[j] & (1<<(i));
      if(t) adj[i] |= 1<<j;
      }
      memset(visited, 0, sizeof(visited));
      }

      void checkmax(int &a, int b) {
      if(b > a) a = b;
      }

      // couting 1's
      int count_bit(int n) {
      int bits = 0;
      while(n) {
      n &= n-1;
      ++bits;
      }
      return bits;
      }



      // bfs
      int search(int status, int N) {
      visited[status] = 1;
      queue<int> Q;
      Q.push(status);
      while(!Q.empty()) {
      int cur = Q.front();
      Q.pop();
      checkmax(maxon, count_bit(cur));
      if(maxon == N) break;
      if(!cur) continue;
      for(int i = 0; i < N; ++i) {
      if( cur & (1<<i)) {
      int dest = cur ^ adj[i];
      if(dest && !visited[dest]) {
      visited[dest] = 1;
      Q.push(dest);
      }
      }
      }
      }
      return maxon;
      }

      int main() {

      //int start = clock();
      while(scanf("%d", &n) != EOF) {
      int N = n*(n+1)/2;
      construct_adj(N); // construct adj matrix

      int cells = 0; // current cell status
      for(int level = 1; level <= n; ++level) {
      int i = level*(level-1)/2; // start point of this level
      for(int j = 0; j < level; ++j) {
      char b; // button status
      scanf(" %c", &b);
      if(b == '.') cells |= 1<< (i + j);
      }
      }
      maxon = -1;
      printf("%d\n", search(cells, N));
      }

      // printf("total time %d", clock()-start);
      return 0;
      }

      Josephus problem

      N人围成一圈,每轮报到k倍数的离开,从下一个人开始重新从1开始报数. 求解最后一个剩下的人最初是第几个。
      问题描述:n个人(编号0~(n-1)),从0开始报数,报到(m-1)的退出,剩下的人继续从0开始报数。求胜利者的编号。

      我们知道第一个人(编号一定是m%n-1) 出列之后,剩下的n-1个人组成了一个新的约瑟夫环(以编号为k=m%n的人开始):
      k k+1 k+2 ... n-2, n-1, 0, 1, 2, ... k-2
      并且从k开始报0。

      现在我们把他们的编号做一下转换:
      k --> 0
      k+1 --> 1
      k+2 --> 2
      ...
      ...
      k-2 --> n-2
      k-1 --> n-1

      变换后就完完全全成为了(n-1)个人报数的子问题,假如我们知道这个子问题的解:例如x是最终的胜利者,那么根据上面这个表把这个x变回去不刚好就是n个人情况的解吗?!!变回去的公式很简单,相信大家都可以推出来:x‘=(x+k)%n

      如何知道(n-1)个人报数的问题的解?对,只要知道(n-2)个人的解就行了。(n-2)个人的解呢?当然是先求(n-3)的情况 ---- 这显然就是一个倒推问题!好了,思路出来了,下面写递推公式:

      令f[i]表示i个人玩游戏报m退出最后胜利者的编号,最后的结果自然是f[n]

      递推公式
      f[1]=0;
      f[i]=(f[i-1]+m)%i; (i>1)

      有了这个公式,我们要做的就是从1-n顺序算出f[i]的数值,最后结果是f[n]。因为实际生活中编号总是从1开始,我们输出f[n]+1

      由于是逐级递推,不需要保存每个f[i],程序也是异常简单:

      DP solution: O(N)
      d[1] = 0
      d[i] = ( d[i-1] + k ) % i
      不光最后一个剩下的符合这个规律,倒数第j个离开的也是这样

      d[i,j] = ( d[i-1, j] + k ) % i

       K = 3, $->  second-last; * -> last
      ---------------------------------
      0 1 2
      $ *

      0 1 2 3
      * $

      0 1 2 3 4
      $ *

      0 1 2 3 4 5
      * $

       

      int d[151];

      // standard joseph problem solution
      // person id is 0-based
      int joseph(int n, int m) {
      d[1] = 0;
      for(int i = 2; i <= n; ++i)
      d[i] = (d[i-1] + m)% i;
      return d[n];
      }

      // another joseph problem extension
      // always kill the 1st guy, then every m-th guy
      // idea: after the 1st guy is killed,
      // the problem (n,m) is changed to a (n-1, m)
      // standard joseph problem
      // the solution's id + 1 is the old problem's id
      int joseph1(int n, int m) {
      int last = joseph(n-1, m);
      return last + 1;
      }

      // the solve subrutine is used to find the exact m
      // to make the 2nd person the last guy to kill
      // in extension joseph1
      int solve(int n) {
      int m = 2;
      while(true) {
      // 0-based index
      int last = joseph1(n, m);
      if (last == 1)
      break;
      ++ m;
      }
      return m;
      }

      Problem:




      • UVa 130 Roman Roulette

      • UVa 133 The Dole Queue



        • just simulate the process

        • using array to keep the status (i tried to use bit vector, but it seems no good for this problem)
                  int live = n;  // live person #
          int p1 = 1, p2 = n; // current position
          while(live > 0) {
          int cnt = 0;
          while(true) { //skip k persons
          if(a[p1]) ++cnt;
          if(cnt == k) break;
          reg(--p2,n)
          }
          cnt = 0;
          while(true) { // skip m persons
          if(a[p2]) ++cnt;
          if(cnt == m) break;
          reg(--p2,n)
          }

          if(p1 == p2) {
          printf("%3d", a[p1]);
          a[p1] = 0; live--;
          } else {
          printf("%3d%3d", a[p1],a[p2]);
          a[p1]=a[p2]=0;
          live -= 2;
          }
          reg(++p1,n);
          reg(--p2,n);

          if(live) printf(",");
          }


      • UVa 305 Joseph



        • We have to check :

        • a1 = (m-1) % n > k

        • a2 = (a1 + m-1) % (n-1) > k

        • a3 = (a2 + m-1) % (n-2) > k

        • ...

        • ak = (a(k-1) + m-1) % (n-k+1) > k
          int gen_data(int k) {
          int n = k << 1; // total number n = 2k
          int m = k+1; // starting to check at k+1
          bool found = false;
          while(!found){
          int ap = 0; // previous (a_i-1)
          int an; // current (a_i)
          for(int i = 0; i < k; ++i) {
          an = (ap+m-1) % (n-i);
          ap = an;
          if(an < k)
          break;
          if(i == k-1 && an >= k)
          found = true;
          }
          if(!found) ++m;
          }
          return m;
          }


      • UVa 402 M*A*S*H

      Monday, September 15, 2008

      Two types of couting change problem (DP)

      1. to minimize the number of coins

      d[i] = min(d[i], d[i - denomination[j]] + 1)
      j
      e.g., coins = {1, 5, 10, 25}, 其实算到100以内足够了,100以上的(n > 100) 用
      (n / 100) x 4 + d[n % 100] 就行了。(整dollar直接4个quarter, 零头查100以内的表)如果还要输出各种coin都用了多少,用个d[i][D] 二维表,第二维直接存各种coin用的数. 上面的例子,

      [0] [1] [2] [3] [4]
      0(min num) 1c 5c 10c 25c
      -------+---------+---+----+----+---|
      |d[1] | 1 | 1| 0| 0| 0|
      |d[2] | 2 | 2| 0| 0| 0|
      ....
      |d[5] | 1 | 0| 1| 0| 0|
      |d[6] | 2 | 1| 1| 0| 0|
      ....
      |d[100]| 4 | 0| 0| 0| 4|


      2. to output the maximum number of ways to count the changes
      Uva judge:

      "The number of ways to change amount A using N kinds of coins equals to:
      • The number of ways to change amount A using all but the first kind of coins, plus
      • The number of ways to change amount A-D using all N kinds of coins, where D is the denomination of the first kind of coin." --- Art of programming contest SE for uva

      Sunday, September 14, 2008

      Lights Out Puzzle

      Problem Intro:
      Related discussion:
      topcoder Problem:
      ZOJ Problem: TBD

      可以把这个问题转化为一个graph 问题,用adjancency matrix(dim = row x col) 可以表示出网格(cell)间的邻接关系。
      每个网格的状态是一个finite field。如果只有on/off,则finite field的运算均mod 2.
      我们在求解中用到的运算有求倒数和求模,所以这里给出了finite field的 invert(), modulate().

      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------
      #define N 10 // maximum vertex #
      #define MAX N*N // maximum edge #

      // --- global variables ---
      int colcount; // integer, number of columns
      int rowcount; // integer, number of rows
      int mod; // integer, number of states of a tile (flip: binary, mod = 2)

      int mat[MAX][MAX]; // integer[i][j]
      int cols[MAX];
      int m; // count of rows of the matrix
      int n; // count of columns of the matrix
      int np; // count of columns of the enlarged matrix
      int r; // minimum rank of the matrix
      int maxr; // maximum rank of the matrix

      // integer[row][col], current states of tiles
      int cells[N][N] = {
      {1,1,1,1,1,1},
      {1,0,2,0,1,0},
      {1,0,0,0,1,0},
      {1,0,0,0,1,0},
      {1,1,1,1,1,0},
      {1,2,2,1,1,1}
      };

      // --- finite field algebra solver
      int modulate(int x) {
      // returns z such that 0 <= z < x ="=">= 0) return x % mod;
      x = (-x) % mod;
      if (x == 0) return 0;
      return mod - x;
      }

      int gcd(int x, int y) { // call when: x >= 0 and y >= 0
      if (y == 0) return x;
      if (x == y) return x;
      if (x > y) x = x % y; // x < y
      while (x > 0) {
      y = y % x; // y < x
      if (y == 0) return x;
      x = x % y; // x < y
      }
      return y;
      }

      int invert(int value) { // call when: 0 <= value < mod
      // returns z such that value * z == 1 (mod mod), or 0 if no such z
      if (value <= 1) return value;
      int seed = gcd(value,mod);
      if (seed != 1) return 0;
      int a = 1, b = 0, x = value; // invar: a * value + b * mod == x
      int c = 0, d = 1, y = mod; // invar: c * value + d * mod == y
      while (x > 1) {
      int tmp = floor(y / x + 0.0);
      y -= x * tmp;
      c -= a * tmp;
      d -= b * tmp;
      tmp = a; a = c; c = tmp;
      tmp = b; b = d; d = tmp;
      tmp = x; x = y; y = tmp;
      }
      return a;
      }
      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------


      在mathWorld中讲解的很清楚,给出初始状态和目的状态,最终找可行路径归结为求解 linear equation.
      求解有三种可能:
      1. 无解
      2. 有唯一解 (矩阵满秩)
      3. 有多个解
      由于这是一个finite field 元素的 linear equation,我们用了相应的gaussian elimination method, O(N^4)

      几个子函数:
      getmat(...), setmat(...) 用来存取adj matrix, 以及作gaussian elimination操作matrix


      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------
      int getmat(int i,int j) { return mat[i][cols[j]]; }
      void setmat(int i, int j, int val) { mat[i][cols[j]] = modulate(val); }
      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------


      initMatrix() 根据问题创建adj matrix,不同问题需要改写这一部分
      cols[] 用来记录消去过程中,列之间的交换


      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------
      void initMatrix() {
      maxr = m < n ? m : n;
      for(int i = 0; i < m; ++i)
      for(int j = 0; j < n; ++j)
      mat[i][j] = 0;
      for (int row = 0; row < rowcount; row++)
      for (int col = 0; col < colcount; col++){
      int i = row * colcount + col;

      mat[i][i] = 1;
      // 4-connected
      if (col > 0) mat[i][i - 1] = 1;
      if (row > 0) mat[i][i - colcount] = 1;
      if (col < colcount - 1) mat[i][i + 1] = 1;
      if (row < rowcount - 1) mat[i][i + colcount] = 1;

      // diagonal connection
      //if (col > 0 && row > 0) mat[i][i-colcount-1] = 1;
      //if (col > 0 && row < rowcount - 1) mat[i][i+colcount-1] = 1;
      //if (col < colcount-1 && row > 0) mat[i][i-colcount+1] = 1;
      //if (col < colcount-1 && row < rowcount-1) mat[i][i+colcount+1] = 1;
      }

      for (int j = 0; j < np; j++) cols[j] = j;
      return;
      }
      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------


      sweep(), sweepStep(), doBasicSweep()构成O(N^4) 高斯消去法
      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------
      void doBasicSweep(int pivoti, int pivotj) {

      if(r != pivoti) {
      for(int j = 0; j < np; j++) {
      int tmp = getmat(r, j);
      setmat(r,j, getmat(pivoti,j));
      setmat(pivoti,j, tmp);
      }
      }

      if(r != pivotj) {
      int tmp = cols[r];
      cols[r] = cols[pivotj];
      cols[pivotj] = tmp;
      }

      for (int i = 0; i < m; i++) {
      if (i != r) {
      int air = getmat(i,r);
      if (air != 0)
      for (int j = r; j < np; j++)
      setmat(i,j, getmat(i,j) - getmat(r,j) * air);
      }
      }
      }

      int sweepStep() {
      int i;
      int j;
      bool finished = true;
      for (j = r; j < n; j++) {
      for (i = r; i < m; i++) {
      int aij = getmat(i,j);
      if (aij != 0) finished = false;
      int inv = invert(aij);
      if (inv != 0) {
      for (int jj = r; jj < np; jj++)
      setmat(i,jj, getmat(i,jj) * inv);
      doBasicSweep(i,j);
      return 1;
      }
      }
      }
      if (finished) { // we have: 0x = b (every matrix element is 0)
      maxr = r; // rank(A) == maxr
      for (j = n; j < np; j++)
      for (i = r; i < m; i++)
      if (getmat(i,j) != 0) {
      //std::cout << "no solution\n";
      return -1; // no solution since b != 0
      }
      return 1; // 0x = 0 has solutions including x = 0
      }
      return -1; // failed in finding a solution
      }

      int sweep() {
      for (r = 0; r < maxr; r++) {
      int status = sweepStep();
      if(status < 0) return -1;
      if(status == 0) return 0;
      if (r == maxr) break;
      }
      return find_min_flip();
      }// ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------


      solve 是整个问题入口,初始化各变量,调用initMatrix生成adj matrix,并填入目标状态
      mat[i][n] (i = 1...n) --- mat矩阵最后一列

      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------
      int solve(int goal) {

      int size = colcount * rowcount;
      m = size;
      n = size;
      np = n + 1;
      initMatrix();
      for (int row = 0; row < rowcount; row++)
      for (int col = 0; col < colcount; col++)
      mat[row * colcount + col][n] = modulate(goal - cells[row][col]);
      return sweep();
      }// ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------

      find_min_flip根据矩阵的秩r,
      如果满秩,直接求出解,否则
      可知有n-r个自由变量,构成mod^(n-r)种组合
      e.g., mod = 2, n-r = 4, 则有2^4 = 16种组合。
      对每一种组合代入原方程组,求解其它变量值,对于所有可能,再统计其中步数最少的输出
      这里用了iterative backtracking, 两重循环,
      外重是当前处理第k个自由变量,
      内重是当前处理的自由变量的第(0...mod-1)种取值可能

      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------
      int find_min_flip() {

      vector sol(n,-1);
      vector best(n, -1);

      if(r == n) {
      int flips = 0;
      for(int j = 0; j < n; ++j) {
      if(j > 0 && j % colcount == 0) cout << endl;
      if(getmat(j,n) > 0) flips += getmat(j,n);
      cout << getmat(j,n) << " ";
      }
      cout << endl;
      return flips;
      }
      else if(r < n) {
      int minflips = n;
      int k = r;
      while (k >= r) {
      while(sol[k] < mod) {
      sol[k] = sol[k] + 1;
      if(k == n-1) {
      for(int j = r-1; j >= 0; --j) {
      int sum = 0;
      for(int jj = j+1; jj < n; ++jj)
      sum += sol[jj] * getmat(j,jj);
      sol[j] = modulate(getmat(j,n) - sum);
      }

      int fl = 0;
      for(int j = 0; j < n; ++j) {
      if(sol[j] == 1) fl++;
      //cout << sol[j] << " ";
      }
      //cout << endl;
      if(fl < minflips) {
      minflips = fl;
      best = sol;
      }
      }
      else k = k + 1;
      }
      sol[k] = -1;
      k = k - 1;
      }
      for(int j = 0; j < n; ++j) {
      if(j > 0 && j % colcount == 0) cout << endl;
      cout << best[j] << " ";
      }
      cout << endl;
      return minflips;
      }
      return -1;
      }// ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------

      test in main(), 6x6 grid, on/off two states, goal state 0 or 1, initial state in cells. 输出的两个矩阵分别对应目标状态为全关和全开,矩阵中为“1”的是应该按的开关。
      // ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------
      int main() {
      mod = 2;
      rowcount = colcount = 6;
      for(int i = 0; i < mod ; ++i){
      solve(i); cout << endl;
      }
      }// ----------------------------------------------------------------------------
      // ----------------------------------------------------------------------------

      outputs:
      1 1 1 1 0 0
      1 0 0 1 0 1
      1 0 0 0 1 1
      1 1 0 0 1 1
      0 0 1 1 0 1
      0 1 1 1 0 0

      0 1 0 0 0 1
      1 1 1 0 1 1
      0 1 1 1 0 0
      0 0 1 1 0 0
      0 1 0 0 1 1
      1 1 0 0 0 1


      For the topCoder problem: LightedPanels. We only need to modify the initMatrix() to support 6-connected adj matrix, then ...

      int LightedPanels::minTouch(vector  board) {

      colcount = board[0].length();
      rowcount = board.size();
      imgcount = 2;

      for (int col = 0; col < colcount; col++)
      for (int row = 0; row < rowcount; row++) {
      cells[row][col] = board[row][col] == '.' ? 0 : 1;
      }

      return solve(1);
      }

      Saturday, September 13, 2008

      O(logn) Fibonacci using Matrix property

      long long divide_conquer_fibonacci(long n) {
      long long i, h, j, k ,t;
      i = h = 1;
      j = k = 0;
      while( n > 0) {
      if(n % 2 == 1) {
      t = j * h;
      j = i * h + j * k + t;
      i = i * k + t;
      }
      t = h * h;
      h = 2 * k * h + t;
      k = k * k + t;
      n = (long) n / 2;
      }
      return j;
      }
      long long data type (-2^63-1 to 2^63-1] can hold until Fib(92).
      1, 1
      2, 1
      3, 2
      4, 3
      5, 5
      6, 8
      7, 13
      8, 21
      9, 34
      10, 55
      11, 89
      12, 144
      13, 233
      14, 377
      15, 610
      16, 987
      17, 1597
      18, 2584
      19, 4181
      20, 6765
      .....
      86, 420196140727489673
      87, 679891637638612258
      88, 1100087778366101931
      89, 1779979416004714189
      90, 2880067194370816120
      91, 4660046610375530309
      92, 7540113804746346429 ----- (19 digits)
      93, -6246583658587674878 ----- (overflow)


      Using the Fast exponentiation to solve linear recurrences, some similar problems follow here.

      timeWatch:

      #include "time.h"
      unsigned int start = clock();
      // .... your code
      std::cout << "Time taken in millisecs: " << clock()-start;

      Monday, June 27, 2005

      cvs from start

      first make a directory to hold the repository, note one repository can include many projects
      cd /home/greeness
      mkdir repo

      add a project named test into that directory
      cd repo
      mkdir test
      cvs init

      checkout the repository to elsewhere, add files into the project and commit the change
      export CVSROOT=/home/greeness/repo
      cd /home/greeness/work
      cvs co test
      cp ../*.cpp test
      cvs add *.cpp
      cvs commit *.cpp

      tunel through ssh to get a copy of the project
      cvs -d :ext:greeness@server.cmu.edu:/afs/blah/repository co modulename

      Saturday, June 18, 2005

      conference reminder

      aamas 06 Hakodate, Japan
      http://www.aamas-conference.org/aamas06-cfp.txt
      * autonomous robots & robot teams

      Important Dates
      ---------------

      Oct 6, 2005: electronic abstract submission deadline
      Oct 9, 2005: electronic paper submission deadline
      Dec 20, 2005: notification

      aaai 2006, July 16–20. Boston

      http://www.aaai.org/Conferences/National/2006/

      ICRA - International Conference on Robotics and Automation
      2006 deadline: 16th Sep, 2005.

      IROS - International Conference on Intelligent Robots and Systems
      2006 deadline: around Feb, 2006

      ECCV 2006
      September 9, 2005 Submission of titles and abstracts
      September 12, 2005 Submission of full papers

      Wednesday, June 08, 2005

      firefox reconfig to speed up

      - input in URL bar: "about:config"
      - find "network.http.pipelining", change its value to "true"
      - find "network.http.pipelining.maxrequests" from 4 to 30
      - add a new integer "nglayout.initialpaint.delay" and set its value to 0

      Thursday, May 19, 2005

      make avi movie from mjpeg logs

      # ./bin/extract_mjpeg ../useful_logs/searchballdemo.4.1.mjpeg
      # mencoder "mf://*.jpg" -mf fps=25 -o output.avi -ovc lavc -lavcopts vcodec=mpeg4
      # mplayer output.avi -vo x11 -loop 0